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CHAPTER  I 
INTRODUCTION 


A.  Purpose  of  Investigation 

Light  propagating  through  the  atmosphere  undergoes  bending  due 
to  variations  in  the  index  of  refraction.  When  this  bending  is 
uniform  one  witnesses  phenomena  such  as  mirages  [1]  and  the  apparent 
flattening  of  the  sun  as  it  approaches  the  horizon.  If  the  index  of 
reflection  is  fluctuating  randomly  the  phase  and  amplitude  of  the 
light  are  perturbed  in  a random  manner.  Phase  fluctuations  give  rise 
to  an  uncertainty  in  the  location  of  a radiating  body  [2].  This 
behavior  is  exemplified  by  the  "dancing"  or  "quivering"  of  the  image 
of  a star  as  Seen  through  a telescope.  Intensity  fluctuations,  or 
scintillations  [3]  as  they  are  commonly  called,  are  typified  by  the 
twinkling  of  stars  or  distant  headlights. 

By  examining  the  temporal  power  spectrum  of  the  scintillations 
one  can  discern  important  properties  of  the  medium  through  which  the 
light  has  propagated.  With  the  ready  availability  of  the  laser  it 
is  therefore  natural  that  we  be  interested  in  the  temporal  scintil- 
lation spectrum  of  a laser  beam. 

A realizable  laser  beam  may  not  always  have  the  perfect  prop- 
erties which  analytical  studies  assume.  The  need  therefore  arises 
to  determine  the  effect  upon  its  scintillation  spectrum  of  a laser 
beam  containing  a certain  amount  of  deterministic  degradation, 
e.g.,  a beam  partially  blocked  by  dirt  specks.  We  are  concerned 
herein  with  the  temporal  scintillation  spectrum  of  just  such  a 
"dirty"  beam. 

Within  this  work  a contemporary  technique  currently  appearing 
in  the  literature  is  employed  in  deriving  a very  general  but  compact 
mathematical  expression  for  the  scintillation  spectrum.  This  ex- 
pression, which  is  restricted  to  the  weak  turbulence  regime,  is 
then  applied  to  the  analysis  of  several  situations  of  great  interest. 
Specifically  the  analysis  is  directed  toward  description  of  the 
effects  of  such  phenomena  as  localized  turbulence,  off-axis  de- 
tectors, and  the  use  of  a laser  beam  which  is  blemished  in  a 
deterministic  sense.  Spectra  obtained  under  these  conditions  will 
be  shown  to  provide  additional  information  about  the  propagation 
medium. 
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This  study  was  motivated  by  a desire  to  provide  theoretical 
support  for  experimentally  obtained  scintillation  spectra.  As 
illustrated  in  Figure  1,  these  spectra  were  digitally  calculated 
from  data  obtained  by  propagating  a laser  beam  over  an  approximately 
one  kilometer  path  between  two  aircraft  flying  abreast.  The  spectra 
thus  generated  displayed  a number  of  pecularities  and  in  the  course 
of  the  effort  to  put  this  behavior  on  a firm  theoretical  basis  a 
number  of  hereto for  unexplored  phenomena  were  encountered. 

I • 
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The  principal  contribution  of  this  effort  lies  in  the  develop- 
ment of  a simplified  formula  for  the  temporal  scintillation  spectrum 
and  its  application  to  several  important  configurations  which  are 
commonly  encountered  in  practical  situations.  The  problem  attacked 
in  this  dissertation  is  therefore  a theoretical  description  of  the 
temporal  scintillation  spectrum  of  a "dirty"  laser  beam. 

The  general  formalism  presented  here  is  also  applicable  to 
arbitrary  transmitter  and  receiver  shapes  such  as  the  beam  transmitted 
from  an  unstable  resonator  or  a Cassegrain  telescope  receiver  or 
transmitter. 
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B.  Perspective 


Whenever  an  electromagnetic  wave  is  propagated  through  the  at- 
mosphere it  is  perturbed  in  a random  manner.  For  radio  frequency 
waves  this  corruption  is  due  to  index  of  refraction  variations  caused 
by  humidity  fluctuations  within  the  atmosphere.  At  optical  fre- 
quencies the  index  of  refraction  variations  are  due  primarily  to 
temperature  fluctuations.  Both  humidity  and  temperature  fluctu- 
ations are  caused  directly  or  indirectly  by  radiant  heating  of  the 
earth  and  convective  heating  of  the  atmosphere  in  conjunction  with 
wind  shears.  The  research  effort  described  herein  in  aimed  primarily 
at  predicting  the  effects  of  turbulent  index  fluctuations  at  optical 
frequencies  but  may  be  generalized  to  longer  wavelengths. 

The  refractive  index  field  is  composed  of  randomly  arrayed 
eddies  of  various  sizes  and  indices  of  refraction.  These  eddies 
have  sizes  I such  that  sl0  < a < Lo  where  t0  and  Lo  are  respectively 
the  inner  and  outer  scales  of  turbulence.  Mechanical  energy  is  in- 
jected by  wind  shears  into  the  large  scale  sizes  which  then  break 
into  smaller  and  smaller  eddies  through  the  phenomena  of  vortex 
stretching[89].  The  energy  thus  cascades  into  smaller  and  smaller 
scale  sizes  until  it  is  dissipated  in  the  form  of  heat  at  scale  sizes 
on  the  order  of  the  inner  scale,  t0. 

If  one  assumes  that  the  turbulent  eddies  are  blown  en  masse 
through  the  wave  transverse  to  the  propagation  path,  then  what  will 
be  observed  at  the  receiver  plane  is  simply  the  moving  diffraction 
patterns  of  eddies  of  various  sizes  and  locations. 

2 

For  short  propagation  paths  (L  < iQ/\  where  A is  the  free-space 
wavelength)  the  incident  wave  undergoes  perturbations  to  its  phase 
only.  As  the  path  length  is  increased,  amplitude  fluctuations  are 
also  introduced.  When  viewed  by  a square-law  detector  (e.g.,  the 
eye  or  a photomultiplier)  the  resultant  power  (irradiance)  fluctu- 
ations are  commonly  referred  to  as  scintillation. 

.By  using  a detector  (photomultiplier)  to  measure  the  variance 
of  the  irradiance  of  a perturbed  wavefront  it  is  noticed  that  to 
a point,  as  the  path  length  or  turbulance  level  increases,  so  does 
the  irradiance  variance.  However  for  very  strong  turbulence  and/or 
long  path  lengths  the  variance  is  seen  to  peak  or  saturate  [4], [5] 
and  to  even  decrease  slightly  [6].  This  effect  has  attracted  much 
attention  partly  because  people  intuitively  expected  the  variance 
to  increase  without  limit. 

Within  the  weak  turbulence  regime  (or  for  moderate  propagation 

lengths)  the  effects  of  index  inhomogeneities  can  be  closely  approxi- 
mated by  the  linear  superposition  of  the  interactions  of  the  coherent 
(unperturbed)  wave  with  turbulent  eddies  of  all  sizes  at  all  locations 
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between  transmitter  and  receiver.  This  is  a single  scattering  model. 

Strong  turbulence  gives  rise  to  multiple  scatterings  which  manifest 
themselves  as  the  aforementioned  saturation  phenomenon.  We  shall  be 
concerned  herein  principally  with  the  weak  turbulence  model.  However 
where  appropriate  modifications  will  be  pointed  out  which  enable 
generalization  to  the  strong  turbulence  domain. 

The  temporal  spectrum  of  the  irradiance  fluctuations  could  easily 
have  been  measured  instead  of  the  variance  of  irradiance.  It  is  just 
this  power  spectrum  in  which  we  are  interested. 


C.  Discussion  of  Previous  Work 

Since  primordial  man  first  turned  his  gaze  toward  the  night-time 
sky  he  has  been  intrigued  by  the  phenomenon  of  scintillation.  (Recall 
the  nursery  rhyme  which  goes  "Twinkle,  twinkle,  little  star  ....")  It 
seems  natural  then,  that  scintillation  be  a topic  of  research  from 
earliest  times.  Further,  it  is  reasonable  that  the  study  of  scintil- 
lation have  its  roots  in  one  of  the  oldest  sciences,  astronomy  [7].  In 
more  recent  times  scintillation  has  been  an  important  device  for  passive 
remote  sensing,  for  example  determining  the  velocity  of  the  inter- 
planetary solar  wind  [8,9],  the  ionospheric  drift  [10,11]  or  terrestrial 
wind  [12].  With  the  advent  of  the  laser,  scintillation  has  become 
important  within  the  realms  of  optical  communication  [13,14]  and  active 
remote  sensing  [15,16,17,18,19]. 

One  of  the  classical  mathematical  tools  with  which  scintillation 
is  predicted,  is  known  as  the  Method  of  Smooth  Perturbations  [20]. 

This  technique,  popularized  by  Tatarski  [20],  hypothesizes  an  iterative 
solution  (for  weak  turbulence)  to  the  scalar  wave  equation  for  the 
natu-al  logarithm  of  the  field.  A great  many  authors  have  employed 
this  method  in  the  study  of  (among  other  phenomena)  scintillation. 

The  simplest  problem  related  to  scintillation  which  has  been 
attacked  with  the  Method  of  Smooth  Perturbations  is  the  calculation  of 
the  log-amplitude  variance.  This  problem  has  been  solved  for  plane 
wave  sources  [21],  spherical  wave  sources  [22,23],  and  the  more  complex 
finite  beam  (viz.  gaussian  beam)  sources  [24,25]. 

Another  topic  routinely  of  interest  is  the  spatial  spectrum  (or  its 
inverse  Fourier  transform,  the  covariance)  of  the  log-amplitude  for 
plane  [26],  spherical  [22],  and  finite  beam  wave  sources  [24]. 

The  next  step  in  complexity  is  the  calculation  of  temporal  spectra 
of,  for  example,  angle  of  arrival  [27,28,29],  phase  [30,31,32]  phase 
difference  [31,33,34]  or  irradiance  [35,36,37,38,39].  The  ter^oral 
scintillation  spectrum  for  plane  waves  was  examined  by  Tatarski[36]. 

An  Interesting  feature  of  his  development  was  that  It  accounted  for  a 
finite  receiver  aperture.  However,  since  the  true  plane  wave  is  merely 
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a mathematical  abstraction,  his  analysis  is  applicable  to  only  a small 
class  of  problems  [40]. 

An  extension  of  Tatarski's  method  to  the  calculation  of  the  spheri- 
cal wave  temporal  scintillation  spectrum  was  provided  by  Clifford  [35], 

In  his  evaluation  he  developed  expressions  for  the  asymptotic  behavior 
of  both  plane  and  spherical  waves.  Solution  of  the  spherical  wave 
problem  represented  a significant  advance  in  the  theory  because  the 
spherical  wave  could  be  closely  approximated  in  practice,  thus  allowing 
experimental  verification. 

The  calculation  of  the  temporal  scintillation  spectrum  of  a 
practically  realizable  source  (the  laser)  was  provided  by  Ishimaru  [39] 
and  later  by  Mironov,  et  al  [37].  Both  authors  presented  numerically 
calculated  spectra  for  a variety  of  laser  beam  configurations. 

The  Method  of  Smooth  Perturbations,  as  it  has  been  used  in  the 
past  is»  howeverr extremely  cumbersome  when  applied  to  an  extended  source. 
Its  use  in  the  analysis  of  each  of  the  aforementioned  problems  necessi- 
tated solving  the  wave  equation  for  a particular  source  field  (plane, 
spherical,  gaussian,  etc.). 

A more  modern  approach  is  to  not  solve  the  wave  equation  for  the 
field  directly,  but  rather  for  the  Green's  function  for  the  particular 
propagation  path  configuration.  The  receiver  plane  field  is  obtained 
merely  by  convolving  the  Green's  function  with  the  incident  field. 

This  extension  of  the  Green's  function  technique  to  propagation  within 
an  inhomogeneous  medium  is  commonly  referred  to  as  the  Extended  Huygens- 
Fresnel  Principle  [41].  It  is  perhaps  paradoxical  that  the  Green's 
function  technique  (for  propagation  within  a homogeneous  medium)  pre- 
dates the  Method  of  Smooth  Perturbations  [42]. 

The  study  with  which  we  are  concerned  employs  elements  of  both 
the  classical  (Method  of  Smooth  Perturbations)  and  the  modern  (Extended 
Huygens-Fresnel  Principle)  techniques.  Virtues  of  each  method  are 
exploited  to  obtain  an  expression  for  the  temporal  scintillation  spec- 
trum of  an  arbitrary  source  field.  This  formula  is  much  simpler  yet 
more- general  than  those  developed  by  Tatarski,  Clifford,  or  Ishimaru. 

In  fact  their  results  are  special  cases  of  our  formula.  This  study 
extends  the  results  of  Ishimaru  and  Mironov  in  that  It  accounts  for  an 
arbitrary  extended  source  and  an  arbitrarily  located  extended  receiver 
aperture. 

A portion  of  our  analysis  deals  with  the  effects  of  extended 
receiver  apertures  and  the  ensuing  "aperture  averaging".  This  topic 
has  also  been  treated  theoretically  by  Fried,  et  al  [43]  and  experi- 
mentally by  Homstad,  et  al  [44]  and  Dumphy  and  Kerr  [45]. 

The  evaluation  of  the  temporal  scintillation  spectrum  for  a 
gaussian  beam  with  an  off-axis  detector  represents  a significant 


contribution  to  the  body  of  knowledge.  Previous  authors  have  always 
dealt  with  axial  detectors. 

Finally,  use  of  our  formula  in  the  modeling  of  a dirty  laser  beam 
and  the, subsequent  evaluation  of  the  scintillation  spectrum  provides 
an  important  description  of  a situation  which  is  universally  encounter 
ed  experimentally,  but  rarely  treated  even  qualitatively.  The  techni- 
ques applied  to  the  modeling  of  a dirty  laser  beam  could  also  be 
applied  to  the  prediction  of  the  scintillation  spectrum  of  a laser 
with  cassegrain  optics. 


D.  Outline  of  Effort 


Chapter  II  provides  the  bulk  of  the  background  plus  the  mathe- 
matical foundation  for  the  remainder  of  the  research  effort.  The 
development  relies  upon  the  extension  of  the  familiar  Huygens-Fresnel 
diffraction  integral  to  the  problem  of  propagation  within  an  inhomo- 
geneous medium.  Use  of  the  Extended  Huygens-Fresnel  Integral  results 
in  an  expression  for  the  temporal  scintillation  spectrum  for  an  arbi- 
trary source  field  in  terms  of  the  generalized  spatio-temporal  second 
order  statistical  moments  of  a spherical  wave.  These  moments  which 
are  functions  of  separations  within  the  transmitter  plane,  the  receiver 
plane,  and  time  are  then  derived.  With  some  easily  justified  approxi- 
mations the  expression  for  the  spectrum  undergoes  tremendous  simplifi- 
cation. The  final  expression,  which  is  very  compact,  has  some  very 
interesting  interpretations.  A description  of  the  final  expression  in 
the  context  of  a phase  grating  model  of  the  atmosphere  completes  this 
chapter. 

The  formula  for  the  scintillation  spectrum  derived  in  Chapter  II 
is  applied  in  Chapter  III  to  the  case  of  a clean  gaussian  beam  source. 
Results  of  the  analysis  are  in  terms  of  analytic  formulae  for  the 
asymptotes  of  the  spectrum  under  a variety  of  conditions.  Such  con- 
ditions include  the  case  in  which  the  detector  is  off-axis  to  the  laser 
beam,  as  illustrated  in  Figure  2,  and  the  case  in  which  the  turbulence 
is  localized.  In  addition  to  asymptotic  results,  numerically  calcu- 
lated spectra  are  presented  for  typical  instances. 

In  Chapter  IV  the  "dirty"  laser  beam  is  studied.  The  dirty  beam 
is  modeled  as  arising  from  the  propagation  of  an  initially  unperturbed 
laser  beam  through  a window  containing  an  imperfection  and  thence  into 
the  atmosphere.  This  situation  is  illustrated  in  Figure  3.  The  window 
imperfection  is  modeled  as  an  axially  located  spot  which  slightly 
shifts  the  phase  of  the  incident  beam.  It  is  demonstrated  that  this 
specific  arrangement  displays  characteristics  of  a much  more  general 
situation.  Within  the  development  the  phase  object  is  approximated 
by  a truncated  series  of  functions  involving  Gauss ian-Hermite  polyno- 
mials. These  functions  are  of  particular  importance  because  they 
closely  approximate  the  modes  of  a laser  and  because  they  are  exact 
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(within  the  paraxial  approximation)  eigenfunctions  of  the  free  space 
wave  operator.  Because  the  polynomial  series  is  truncated,  a procedure 
is  demonstrated  which  optimizes  the  approximation. 

The  expression  resulting  from  the  modeling  of  a dirty  laser  beam 
is  then  numerically  evaluated  for  a variety  of  situations.  As  a result 
of  the  analysis,  the  scintillation  spectrum  is  shown  to  exhibit  very 
interesting  and  informative  low  frequency  behavior. 

The  summary  and  conclusions  are  contained  in  Chapter  V.  Sug- 
gestions are  provided  for  future  analytical  as  well  as  experimental 
efforts. 


9 


CHAPTER  II 


DERIVATION  OF  FUNDAMENTAL  FORMULA 


A.  Introduction 

This  chapter  is  devoted  to  the  derivation  of  an  expression  for 
the  temporal  scintillation  spectrum  of  an  arbitrary  incident  field. 

The  final  formula  obtained  is  the  foundation  on  which  the  remainder 
of  the  analysis  rests. 

Within  the  derivation  use  is  made  of  the  Extended  Huygens-Fresnel 
principle  and  the  weak  turbulence  approximation.  The  final  result  is 
a very  simple  and  physically  interpretable  formula  for  the  temporal 
scintillation  spectrum  of  an  arbitrary  source  field.  This  formula, 
which  requires  knowledge  of  the  free-space  receiver-plane  field  of 
an  unspecified  extended  source,  also  accounts  for  an  extended  receiver. 

In  Section  B we  employ  the  Extended  Huygens-Fresnel  Integral  to 
develop  an  expression  for  the  scintillation  spectrum  in  terms  of  the 
second  order  spatio-temporal  statistical  moments  of  spherical  waves. 
These  moments,  (log-amplitude  covariances,  phase-1 og-ampl itude  cross 
covariances,  and  wave  structure  functions)  which  are  functions  of 
separation  in  the  transmitter  plane,  separation  in  the  receiver  plane, 
and  separation  in  time,  are  derived  in  Section  C by  use  of  the  method 
of  smooth  perturbations. 

The  results  of  the  two  previous  sections  are  combined  in  Section  D. 
Use  is  then  made  (and  justified)  of  the  weak  turbulence  approximation. 

By  rearranging  terms  of  the  resulting  formula  we  obtain  a series  of 
diffraction  integrals  which  are  then  performed  symbolically  to  yield 
the  desired  expression  for  the  scintillation  spectrum. 

In  Section  E the  various  components  of  the  final  formula  are 
interpreted  in  the  context  of  a phase  screen  model  of  the  atmospheric 
turbulence. 

Finally  Section  F consists  of  a discussion  and  summary  of  the 
chapter. 


B . Derivation  of  Genera  1 Expression 
for  Scintillation  Spectrum 

This  section  is  devoted  to  developing  an  expression  for  the  tem- 
poral scintillation  spectrum  in  terms  of  the  second  order  statistical 
moments  of  a spherical  wave.  The  derivation  is  based  upon  use  of  the 
Extended  Huygens-Fresnel  integral  and  employs  standard  mathematical 
procedures  currently  employed  in  the  literature. 

We  are  interested  in  the  situation  in  which  an  arbitrary  unper- 
turbed source  field  is  incident  (at  z=0)  upon  a region  (z>0)  contain- 
ing randomly  varying  index  of  refraction  inhomogeneities.  The  field 
then  propagates  through  the  inhomogeneous  medium  a distance  L to  a 
receiver  plane  in  which  is  located  an  arbitrary  detector.  To  deter- 
mine the  receiver  plane  fields  we  propose  to  solve  the  scalar  wave 
equation 


[v2  + k2n2(R)]E(R)  = 0 


(1) 


where  the  index  of  refraction  n is  a random  function  of  space  and  k 
is  the  free-space  wavenumber  (k  - 2tt/a).  By  defining  a generalized 
Green's  function  such  that 


[v2  + k2n2(R)]G(R,R')  = - 4*6(|R-R|) 


(2) 


it  is  easily  shown  (see  Appendix  A) [41]  that  the  receiver  t lane  field 
may  be  expressed  approximately  as 


E<R'>  -T7 


r G(r,R')E(r)dr 


(3) 


where  the  surface  integration  is  carried  out  over  the  plane  z=0  and 
L i-s  the  distance  between  transmitter  and  receiver  planes.  In 
Equation  (3)  we  have  denoted  the  three-dimensional  vector  by  an  upper 
case  letter  and  the  two  dimensional  vector  by  a lower  case  letter. 

We  shall  attempt  to  retain  this  convention  throughout  the  remainder 
of  this  work.  Using  the  explicit  expression  for  the  generalized 
Green's  function  (also  discussed  in  Appendix  A)  Equation  (3)  may  be 
wri tten 


E(ri,t)  = Sr  | drl  e 


_ ikl(r1,r|)  + ^(r^rjit)  _ 


E(ri) 


(4) 


where  L is  the  geometrical  (free  s pace )_di stance  between  the  point  r 
in  the  transmitter  plane  and  the  point  r'  in  the  receiver  plane  and 
ip  is  a complex  phase  perturbation  to  the  field  of  a spherical  wave. 
The  variables  of  interest  are  illustrated  in  Figure  4. 


Figure  4--Variables  appearing  in  Huygens-Fresnel 
formula. 


The  preceding  discussion  aimed  at  solution  of  the  wave  equation 
within  an  inhomogeneous  medium  encompasses  the  method  commonly  referred 
to  as  the  Extended  Huygens-Fresnel  technique.  This  topic  is  discussed 
more  fully  in  Appendix  A. 

Using  Equation  (4)  we  express  the  irradiance  at  r'  as 


Ifrjltj) 


|E(rj ;t-| ) |^  = jdrl  ldr2 


ik[L(r1 ,rj)-L(r2,rj)] 
e 


*(n;r1;ti) + ♦*(r2;n,;ti) 

e E(r])E*(r2) 
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Note  that  we  have  indicated  the  covariance  of  irradiance_to_be  simul- 
taneously homogeneous  in  the  receiver  plane  separation,  r-J -rA > and 
stationary  in  time,  t.  This  assumption  will  be  justified  later  in 
this  chapter. 

Use  will  now  be  made  of  Equation  (5)  to  calculate  the  irradiance  cor- 
relation and  mean  irradiance  as  required  in  Equation  (11-b). 

The  mean  irradiance  from  Equation  (2)  is  simply 

ik[L(r1,r])  - L(r2, rj)] 


<I(ri;V>  -(stT  fdri  |dr2 

x E(r])E*(r2)  <e^> 
where  we  have  defined 

c = <l»(n  ;rj  ;t-|)  + **(r2;r|  ;t]) 


(12) 


(13-a) 


If  we  define  \p  = x + iS,  where  S is  the  phase  and  x is  commonly  called 
the  log-amplitude  then  z,  is  given  by 


£ = x(r>| sr-J ;t-j ) + iS(r-| ; r-J ; t-, ) + x(r9;rj;ti)  - iS(r?;rj  ;t-, ) . 


1 * ' 1 ’ **1 


2 ’ ' 1 ’ L1 


2,r  1 * L1 


(13-b) 


The  notation  in  the  derivation  to  follow  is  greatly  simplified  by 
adopting  the  convention 


c = (Xm  + *211)  + i(sm  " s2i l ) 


(13-c) 


where  the  subscripts  refer  to  respectively  a point,  r-j  or  r2,  in  thp 
transmitter  plane,  a point,  r-j  or  r£,  in  the  receiver  plane,  and  a 
time,  t-j  or  t2. 

Under  the  assumption  that  the  complex  phase  perturbation  is  normally 
distributed  [46]  we  have  [47] 


5 + U 

<e^>  e 4 c 


(14) 
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where  op  and  y are  respectively  the  variance  and  mean  of  c-  The 
mean  isc 


= 2<x>  . ( 1 5-a ) 


It  is  easily  demonstrated  by  conservation  of  energy  arguments  [48] 
that  the  mean  must  be  the  negative  of  the  variance; 


= 2<x>  = -2CX(0;0;0) 


(15-b) 


For  an  infinite  plane  wave,  relation  (15-b)  is  identically  true  be- 
cause the  average  irradiance  at  a point  in  the  receiver  plane  must  be 
independent  of  turbulence  strength  (energy  must  be  conserved).  For 
a spherical  wave  propagating  through  weak  turbulence,  energy  is  not 
actually  diffracted  out  of  the  beam  but  merely  instantaneously  redis- 
tributed. Therefore,  for  weak  turbulence,  the  mean  irradiance  of  a 
spherical  wave  should  also  be  independent  of  turbulence  strength. 

The  variance  of  s is  given  by 


p O 

= <C(xm  - <x>  + X2H  - <x>)  + i(s-in  - s211)r> 

= 2 Cx(0;0;0)  + 2 Cx(?r?2;0;0)  (16-a) 

- i4<x><(S111  - S211)>  + i2<(Xlll  + x2n)(s111  - s211)> 

" DS  (r-|-r2;0;0)  (16-b) 


where  is  the  phase  structure  function  defined  as 


Ds(p;p';t)  = 2[CS(0;0;0)  - Cs(p;p';t)] 


(16-c) 


and  we  have  assumed,  simultaneous  homogeneity  in  both_the  transmitter 
plane  separation,  p,  and  receiver  plane  separation,  p1,  and  station- 
arity  in  the  time  separation,  r.  This  assumption  of  simultaneous 
homogeneity  will  be  justified  in  a later  section. 

The  imaginary  terms  on  the  right  hand  side  of  Equation  (16-b)  are 
identically  zero.  This  can  be  verified  by  interchanging  the  points  r^ 
and  r2  in  the  transmitter  plane  thus  introducing  a sign  change  in  these 


'Wjp 


terms.  However  by  isotropy,  this  interchange  should  produce  no  sign 
change.  Therefore  the  imaginary  terms  are  zero  and  Equations  (14), 
(15-b)  and  (16-b)  yield 

-%Ds(ri-r2;0;0)  - ^(rj-rg ;0;0) 

<e?>  = e (17-a) 


or 


<e 


c>  = 


-%0w(r]-r2;0;0) 


07-b) 


where  Dy  is  the  "wave"  structure  function  [48]  defined  simply  as  the 
sum  of  the  phase  and  log-amplitude  structure  functions. 

The  mean  irradiance  term  in  Equation  (11-b)  is  then  finally 
<l(rj;t,)><I(r$;t2)>  . |df,  |df2  JdF3  JdF4 

x E(n)E*(r2)E(r3)E*(F4)  .rp-£(F2.Fi 


Dy  ( ri  ” r2 » 0 ; 0 ) + 0w(r3-r4;0;0)] 


(18) 


Similarly  the  correlation  of  irradiance  in  Equation  (11-b)  is 
given  by 


:I(rj;tj)I(r£st2)>  - (5^  [dr,  |dr2  Jdr3  jdr4 
:E(F1)E*(F2)E(F3)E.(F4)6,kU<?’>fi)-t(F2^)H(F3'?2)-L'F4-f2>] 


x <e^> 


(19-a) 


where 


£ = *(r1;rj;t1)+**(r2;r];t1)+*(r3;r^;t2)+i|»*(r4,r^;t2)  (19-b) 
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By  the  previously  adopted  notation  we  have 


£ = (Xm  + x211  + x322  + x422)  + KSm  - S2H  + S322  - S422)  • 

(201 

Reapplication  of  Equation  (14)  to  calculate  the  mean  of  exp(sl- 
requires  the  calculation  of  the  mean,  and  the  variance,  o|,  of  £. 

Again,  conservation  of  energy  [48]  requires  that  the  mean  be  equal  to 
the  negative  of  the  variance; 


= 4<x>  = -4  Cx(0;0;0) 


(21) 


The  variance  is 


ac  = fi  + f2  + f3  + U 


(22-a) 


where 


fl  = <(xill  “ <x>  + X211  ■ <X>  + X322  " <X>  + X422  " <X>)  > 


f2  = -is  <x><(sm  - s2n  + S322  - s422)> 


(22-b) 

(22-c) 


^3  = 12  <(xm+  x2n  + X322  + X422)^slll  • s211  + S322  * S422^ 

(22-d) 


and- 

f4  = " <^slll  “ s211  + s322  “ s422^2>  (22-e) 

The  term  f-|  is  easily  shown  to  be 
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f]  = 4 Cx(0;0;0)  + 2 Cx(r]-r2;0;0)  + 2 Cx(r3-r4;0;0) 

+ 2 Cx(r1-r3;r1-r^;t1-t2)  + 2 Cx(r|-r4;r|-r$;ti-t2) 

+ 2 Cx(r2-r3;rj-r^;t1-t2)  + 2 Cx(r2-r4;rj-r^;t1-t2)  . (23) 

Tenn  f2  is  identically  zero  because  of  stationarity  and  term  f3  with 
isotropy  arguments  similar  to  those  employed  previously  is  given  by 


f3  = i 2[<(S322  “ S422) (xi l 1 + x21l)> 

+ <(Sir|  - S2i 1 )(x322  + x422)>] 


(24-a) 


f3  = i4[CxS(rrr3;r^r^;trt2)  - CxS(r2-r4;rj-r ' ;trt2)] 


(24-b) 


where  Cx3  is  the  cross-covariance  between  phase  and  log-amplitude. 
Use  of  the  identity 

(a-b+c-d)2  = (a-b)2  + (c-d)2  + (a-d)2 


+ (b-c)^  - (a-c)2  - (b-d ) 2 


gives  for  the  remaining  term  in  Equation  (22-a) 


f4  - -[<(Sm  - S2n)2>  + <(S322  - S422)2>  + <(5111  - S422)2> 


+ <(S2n  - S322)2>  - <(S]n  - S322)2>  - <($21 1 " s422)2>] 


(26-a) 


= "C^S(ri"r2’0’0)  + Ds(r3-r4;0;0)  + Ds(rrr4;rrr2;trV 
+ Ds(72-73;7]-7';trt2)  - Ds(7r73;7j-7' ;trt2) 

- Ds ( ^2"^4 "^2 ; tr ^ (26-b) 

Finally,  combining  Equations  (19)  through  (26)  gives 

<e^>  = exp  |CDw(71-72;0;0)  + M^r^050) 

+ Dw(7r74;7j-7^trt2)  + Dw(72-73;7j-7^;t1-t2) 

- Dw(r-|-r3;r]-r2;t-|-t2)  - Dw(r2-r4;r-j-r2;t-|-t2)] 

+ 2[Cx(r1-73;rj-r^;t1-t2)  + Cx(72-74;7|-7£;ti-t2) 

+ i CxS(71-73;7j-72;t-,-t2)  - i Cxs(72-74;7i-7^;t1-t2)]  \ . 

J(27) 

The  optical  distance  term  in  Equations  (18)  and  (19-a)  is 

ik L = ik[L(71,7j)  - L(72,7j)  + l(F3, r£)  - L(74,7£)]  • (28-a) 


Referring  to  Figure  5 it  is  seen  that 


iki.  = ik[(L2+|r-|-r-j  I2)*5  - (L2+|r2-rj  I2)1® 
* (L2t|F3-r^|2)1*  - (l2+|F4-^|2rt 


(28-b) 


Invoking  the  paraxial  ray  approximation  (|r^-rj|  « l 4M,j)  we  obtain 
iki-  = Cl 1 2 - |r2-ri|2  + |7rri|2  - |74-7i|2] 


TRANSMITTER 

PLANE 


RECEIVER 

PLANE 


k 


Figure  5— Statistical  moment  vector  arguments. 


We  now  obtain  by  combining  Equations  (10),  (11),  (18),  (19-a), 
(27),  and  (28-c)  the  basic  expression  for  the  temporal  Irradlance 

I spectrum; 

sp(“>  ' (kl)4  L d'  e'tm 

x jdr-i  |dr2  jd?3  |dr4  E(F1  )E*(F2)E(F3)E*( F4) 

elr nvii2  • ivii2  + ivii2  - ivii2] 

Ix  H(r|,r2,r3,r4;rj,r£;T)  (29-a) 

where 

K 

I 
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H = exp  |-  l[fy(F,-F2;0;0)  + Dy (F3-F4 ;0;0) 

+ + Dw(F2-F3;Fj-F^;T) 

- Dw(rrF3;Fj-F^;T)  - Dw(F2-F4;F]-F^jt)] 

+ 2 [ Cx  ( Fj  -F3  ;Fj  -F2 ; t ) + Cx(F2-F4;Fj-F^;T) 

+ i CxS (Fj -F3 ;F| -F2 ; t ) - 1 CxS(F2-F4;F]-F^;t)]| 

- exp  .["-^{Dy^-F^OiO)  + Dw(F3-F4;0;0)]"l  . (29-b) 


At  this  point  it  may  appear  that  the  spectrum  is  not  real  as  it 
should  be.  This  fear  may  be  removed  by  manipulating  the  dunmyvarl-  _ 
ables  in^  Eqs.  (29).  Specifically  by  interchanging  r^  and  r4,  r2  and  r3; 
rl  and  r2;  and  setting  t=-t  the  right  hand  side  of  Eq.  (29)  appears  to 
have  been  conjugated.  However  since  the  left  hand  side  is  unchanged* 
the  spectrum  S is  equal  to  its  complex  conjugate,  i.e. • it  is  real. 

Equations  (29)  are  the  final  result  of  this  section.  The  particu- 
lar form  of  the  transfer  function  H is  important  because  at  least  one 
group  of  workers  [49]  has  found  that  for  very  strong  turbulence,  the 
behavior  of  the  H function  is  dominated  by  that  of  the  structure 
functions.  That  is,  the  covariances  in  Equation  (29-b)  may  be  neglected. 
Further,  under  certain  circumstances  (on-axis  point  detector),  approxi- 
mate closed-form  expressions  for  the  wave  structure  functions  may  be 
found.  For  this  instance  it  would  be  possible  to  calculate  the  temporal 
scintillation  spectrum  under  saturation  conditions.  This  topic  seems  a 
fertile  area  for  further  research. 

We  have  in  this  section  derived  an  expression,  in  terms  of  the 
second  order  moments  of  a spherical  wave,  for  the  temporal  scintilla- 
tion spectrum  of  an  arbitrary  field  source.  Use  was  made  of  the  Extended 
Huygens-Fresnel  Integral  which  expresses  the  receiver  plane  field  as 
the  convolution  of  the  transmitter  plane  field  with  a Green's  function 
which  has  been  generalized  to  propagation  within  an  inhomogeneous  med- 
ium. With  standard  mathematical  procedures  and  assumptions  treated 
extensively  in  the  literature  we  arrive  at  the  desired  formula  expressed 
in  Equations  (29). 
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Finally,  approximations  were  discussed  which  may  enable  the  use 
of  the  final  expression  to  calculate  the  scintillation  spectrum  within 
the  saturation  regime. 

The  next  section  is  devoted  to  a derivation  of  the  statistical 
moments  required  by  Equation  (29). 


C.  Derivation  of  Statistical  Moments 

We  now  face  the  task  of  deriving  expressions  for  the  statistical 
functions  required  by  Equation  (29-b)  of  the  previous  section.  Use 
will  be  made  of  the  method  of  smooth  perturbations  [2]  which  is  some- 
times referred  to  as  the  Rytov  method  [50,51].  Since  this  method  is 
well  known,  only  the  salient  features  of  the  derivations  are  presented 
within  this  section.  The  expressions  for  the  statistical  moments 
which  are  developed  within  this  section  represent  generalizations  of 
expressions  currently  available  in  the  literature. 

The  method  of  smooth  perturbations  relies  upon  the  expression  of 
the  field  as  [2] 

l A, 

E = e n=0  (30-a) 

where 

Y =y<[ n(R)  - <n(R)>]2> 

is  the  R.M.S.  variation  in  the  index  of  refraction. 

Typically  y is  on  the  order  of  10"®  [41].  By  inserting  E from  Equation 
(30-a)  in  the  scalar  wave  equation  and  equating  like  powers  of  y one 
arrives  at  an  infinite  series  of  linear  constant  coefficient  (Recall 
that  the  wave  equation  within  an  inhomogeneous  medium  has  the  variable 
coefficient  n2(R).)  differential  equations  which  are  easily  solved  by 
the  Green's  function  technique.  The  solutions  of  these  differential 
equations  are  commonly  denoted  as  <|>n  and  are  given  by 


ll i = vn4' 

vn  T n 


(31) 


Obviously  this  results  in  the  field  being  given  by 
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00 

l *n 

£ = e n=0 

(32) 

Now  if  we  define 

E = a eiS  = e*nA+iS 

(33) 

and  denote 

+n  = *nAn  + iSn 

(34) 

we  see  that 

oo 

*n(A/A0)+i(S-S0)  n=l  *" 

e = e 

(35) 

The  log -amplitude  and  phase  perturbations  which 
respectively  as 

are  denoted 

x = *n(A/A0) 

(36-a) 

and 

tO 

II 

to 

to 

o 

(36-b) 

are  then  to  first  order  smallness  of  y; 

x = Re{i|/]} 

(37-a) 

and 

S-j  = Im{ij>i } 

(37-b) 

Tatarski  [52]  gives  for  the  first  perturbation 

to  the  complex  phase 

f eiklR~RI 1 E(R) 

| |R-R'|  E(R’) 


nT (R)dR 


(38-a) 
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where 


n^R)  = n(R)  - <n(R)>  (38-b) 

and  the  integration  is  carried  out  over  the  scattering  volume  denoted 
by  V.  Figure  6 illustrates  the  variables  used  in  Equation  (38-a). 


4»-|  ( R‘  ;t)  = (r^  ;t)  = j dR  n^Rjt) 


ik[|R'-R|+|R-ri |-|R,-ri H 

■ K ~ri ' e 

x — — — (40) 

|R-r.|  | R 1 -R | 

where  we  have  indicated  explicitly  the  time  dependent  behavior  of  n-j. 
This  time  dependence  arises  because  we  assume  the  wind  to  be  blowing 
the  index  inhomogeneities  across  the  propagation  path. 

Utilization  of  the  paraxial  approximations 


|R-ril  = z + 


(x-xi)2+(y-y.)2 


. (xj-x.)2+(yj-yi)2 

I R,_ri I = L ♦ 1 2L 


(xj-x)2+(yi-y)2 
lR'-Rl  ■ (L-z)  ♦ ~ L2TL-z) 


allows  Equation  (40)  to  be  written 


= r,  £dz  (iXpt)  Dxl..dy  "i(x-y>Zit) 


(xj-x)2+(y^-y)2  (x-xi)2+(y-y1)2  (x^-Xj)2+(y^yj)2 


x exp<ik 


In  Equation  (42)  the  magnification  terms  have  been  approximated  by  the 
first  term  of  their  Taylor's  series  expansions: 
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V 


IR’-r. 


|R-F.||R'-R|  z(L-i) 

To  examine  Equation  (42)  in  the  spatial  spectral  domain  we  make 
use  of  the  Fourier-Stieltjes  expansions  [53] 


r°}  i(<;x+Kly) 

n-j (x,y,z;t)  = IJ  e dv(ie^  ,tc‘  ,z;t) 


(44-a) 


QD 

h<Fi;t)  * If 


d*( <i »<2 » L » t) 


so  that 


(44-b) 


II  ■ s J*  (&)Ildx  I..dy 


x e 


i(<^x+K^y) 


dv(<j,K^,z;t) 


x exp 


Rx]-x)2+(y^-y)2  (x-x.)2+(y-y.)2  (x-j-x. )2+(y^-yi )‘ 


2(L-z) 


In  Appendix  B it  is  demonstrated  that  inversion  of  Equation  (45) 
yields  for  the  spatial  spectrum  of  the  first  perturbation  to  the  com- 
plex phase; 


where  v is  the  spatial  spectrum  of  the  index  of  refraction  fluctuations. 
Since  the  expansion  of  Equation  (44-b)  takes  place  in  the  receiver 
plane,  the  variables  k-j  and  correspond  to  spatial  scales  within  the 
receiver  plane. 

Because,  by  definition,  the  log-amplitude,  x,  is  real,  it  is 
easily  demonstrated  [54]  that  (to  the  first  order  smallness  in  y) 

F i 

X - <X>  ,K2>L;t)  + d**(-<1  ,ric2,L;t)] 

= da(K1(K2.L;t)  , (47) 

where  F indicates  the  two  dimensional  spatial  Fourier-Stieltjes  trans- 
form. 

Using  the_fact  that  the  index  of  refraction  field  is  real  (thus 
requiring  dv(-ic,z;t)  = dv*(«,z;t)) , Equations  (46)  and  (47)  give 


Similarly,  to  the  first  approximation,  the  phase  fluctuation  is 
given  by 
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which  is  easily  shown  to  give 


rL  (\  xiikiilr..K 

C1  tK2,L;t)  = kj  dzdv  *1 , - <2,z;tje  z 


fklkil  2' 


cos  2kF^K 


As  indicated  in  Equations  (47)  and  (49)  we  have  made  use  of  the 
Fourier-Stieltjes  expansions  of  the  log-amplitude  and  phase; 


x(ry,rj;t)  - <x> 


i (<1 X1 +*2yi ) 


da • (<i  ,*2,L;t) 


(51a) 


- - f(°°  i(Kixi+K2yi) 

^1  (ri ; r-J  It)  = e d4>i(<1  ,tc2,L;t) 

• • — oo 


(51b) 


where  the  subscript  i on  the  spectra  denotes  the  source  point  r . . 
Now  the  time  delayed  log-amplitude  and  phase  are  given  by 


X(rj;r2;t+T)  - <x> 


i^x^+K^) 


daj(<.|  ,<2;t+T)  (52a) 


— — f f°°  i(K1x2+K2^2^ 

S1 (rj ;r2;t+T)  = JJ  e d^j (<1 ,<2 ;t+x) 


(52b) 


Introduction  of  Taylor's  frozen-turbulence  hypothesis  [55]  which  says 
eddies  do  not  change  shape  in  the  time  it  takes  them  to  blow  through 
the  propagation  path. 


0|(r,z;t+T)  - n^r-VT.z-.t), 


where  r = (x,y)  and  v is  the  component  of  the  wind  velocity  transverse 
to  the  propagation  path,  results  in  the  spectrum  of  Equation  (52a) 
being  given  by 
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-"‘“—‘“."■-IT 


daj(Kl  »K2»*-;t+T)  = k dzdv^Ki»  \ <2>z>t) 


(54) 


i 


The  expression  for  the  phase  spectrum  is  identical  except  the  sin 
function  is  replaced  by  a cos  function. 

Since  subsequent  derivations  for  the  desired  statistical  quanti- 
ties are  all  similar,  we  will  treat  explicitly  only  the  log-amplitude 
covariance  and  state  without  proof  expressions  for  the  phase  and 
phase-log-amplitude  covariances  and  the  wave  structure  functions. 

The  log-amplitude  covariance  defined  as 

= (55) 

<[x(ri;r^;t1)  - <x>][x(rj;r£;t2)  -<x>]> 


is,  from  Equations  (51a)  and  (52a), 


cx(rrrjiri'r2;T)  = 


(56) 


i < - r’  cp  -i < ' - rA 

e j e ^^a^^  ,<2,L;t)dat(K1  ,<2,L;t+T)>  . 


With  Equations  (48)  and  (54)  the  covariance  is 
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r 00 


r r iK>ri 

cx(rrrjiri"r2;T>  * JJ_„ e 


x k‘ 


L rL 

dzl 

„ 1 


dZp  sin 


L(L-Zi ) 2 


2kz, 


sin 


LO^)  2 

2kz„ 


(57) 


x exp 


i (L-z-j ) 


i(L-z2) 
z2  " 


rj' 


x<dv(z1  V z7  K2‘zl;t)dv*(z2  K\*  z2K2,Z2;t)>- 


The  assumption  of  wide-sense  stationarity  of  the  turbulence 
provides  the  (statistical)  orthogonality  of  the  index  spectra  ex- 
pressed as  [56] 


where  Fn  is  the  two-dimensional  power  spectrum  of  the  index  of 
refraction  field. 

Combining  Equations  (57)  and  (58),  making  the  changes  of 
variables 


and 


and  performing  the  integrations  over  ui  and  U2  results  in 


(59) 
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cx<Vrj;rrrrT>  - k?lf  . de  £ dzi£  dz; 


z,  (L-z^ ) pi 

x Fn(<i  ,k2,  Izj-Zgl ) sin  — 2j^[ — < sin 


z2^~z2^  2 


1 k • * t,  *(il  rj  - £ >|)+  r,  - -^1  rjjj  , 


x exp  i 


where  we  have  dropped  the  primes  on  the  spatial  frequencies. 

Recall  that  in  Equation  (46)  the  spatial  frequencies  correspond 
to  scale  sizes  at  the  receiver.  However  the  changes  of  variables  in 
Equations  (59)  produce  in  Equation  (60),  spatial  frequencies  < cor- 
responding to  spatial  inhomogenieties  within  the  medium.  This  is  an 
important  distinction  for  the  subsequent  manipulation  of  Equation 
(60). 

We  shall  now  proceed  to  simplify  the  expression  for  the  log- 
amplitude  covariance.  The  method  will  make  use  of  a sum  and  dif- 
ference change  of  (path  integration)  variables  and  several  order- 
of-magnitude  arguments. 

The  kappa  integrand  of  Equation  (60)  is  of  the  form 


fL  fL 

dz,  dZo  G(k,z, )Zn) > 

Jn  1 J n c 1 c 


which  upon  changing  to  sum  and  difference  coordinates 


zl+z2 

z = — — » z1  = z + P/2 


p = z1-z2  , z2  = z - p/2 


yields 
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fL  fL  fL/2  (Iz 

j dz-J  dz2  G(k,z1,z2)  = J dzj  dP  G(*,z+p/2, z-p/2) 


f2(L-z) 


(63a) 


+ [ dz  [ dp  G(k,z+p/2,  z-p/2), 

JL/2  J-2(L-z) 


where  the  G function  is  given  by 


G(k,z+p/2, z-p/2)  = Fn(K1,K2,|p| )e 


1 TIC  • V 


(63b) 


x sin 


(z+p/2)(L-z-p/2)  2 


2kL 


I 


x expj^-  • [ (z+p/2)r-j- (z-p/2 )r2+(L-z-p/2)ri-(L-z+p/2)r-] 


-h 


Using  the  trigonometric  identity  for  the  product  of  two  cosines  and 
factoring  the  exponents  gives  for  the  G function 


G(k,z+p/2, z-p/2)  = Fn(ic,  ,k„,|p|  ) 


nv  1 *2 

x exp^lK.|f  (ri-fjj+i  I - fl^-r,)  + xv 

l-LlMi-frl  JL 


<ri-rp+(i  - r)(rrrj: 


(64) 


x e 


JJI—  - I K2  - 

1 k . 


>1 

x ^cos|r  'gkL-'-  K~  "cos 


For  an  isotropic  field,  fluctuations  separated  a distance  p 
are  correlated  only  by  index  of  refraction  inhomogenieties  of  size 
t such  that  [57] 


P £ * 


(65) 


Since  the  dimension  t corresponds  to  a spatial  frequency 

t = 2*/k, 


32 


the  two  dimensional  index  of  refraction  spectrum  Fn  is  markedly 
different  from  zero  only  for 


p £ 2if/K 

(66a) 

tcp  S 2tt. 

(66b) 

Within  the  region  of  significance  of  F , the  argument  of  the 
first  cos  function  in  Equation  (64)  is  n 


p ^2k^Z"^'  ^ ^(1~2z/L)k  <_  * 


A/A 


(67) 


where  *0  is  the  inner  scale  of  turbulence  (typically  V''10~3-10"2m). 
Because  the  wavelength  at  optical  frequencies  (xv|0-6m)  is  much  smaller 
than  the  inner  scale,  we  may  let 


cos  p 


2]%  1. 

Also  within  the  aforementioned  region 


(68) 


k ' [(rj+r^-fryrj)]  if  I (rj+r^)-(ri+rJ.)  I . (69) 

If  we  restrict  our  attention  to  lateral  distances  from  the  propa- 
gation axis  much  smaller  than  the  range  - the  paraxial  approximation 
has  already  required  this  - we  can  allow 


2T  K • [(r-j+rAMr.+r.)] 
e c J * 1 


(70) 


With  the  preceding  approximations  we  now  obtain  for  the  in- 
tegrand of  Equation  (63a), 


G(k,z+p/2,  z-p/2)  = Fn(<1,K2,|p|)sin2 


x exp 
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by 


Finally,  the  log-amplitude  covariance  in  Equation  (60)  is  given 


" f 1 /2  2z 

Cx ( r*i -rvj ; rj - ; T ) = 2k2L2  jj  dic  j{0  dz  { Qd|5  GX(K »z+pL/2,z-pL/2) 


r 1 r2(l-z) 

+ dz  dp  6 (K,z+pL/2,z-pL/2) 

J 1/2  Jo  x 


(72-a) 


where 


G * F (k  k oL)  sin2  f211'21-1-  *2  - Lt«p)H 
fax  V i’  2,p  ' sin  2k  8k  J 

x exp{i7  * [z(rj-r2)  + (l-z)(F1-FJ- ) + tv]}  , (72-b) 

and  we  have  replaced  the  z variable  of  integration  by  z/L. 

It  is  easily  demonstrated,  using  arguments  identical  to  the  pre- 
ceding, that  the  phase  covariance  and  phase-log-amplitude  cross- 
covariance are  given  respectively  by 

00  ( 1/2  2z 

CsCFi-F^Fj-F^T)  = 2k2L2  jjdF^J  dz  | dp  Gs(ic,z+pL/2,z-pL/2) 


r1  f2(l-z) 

^ ^dz  dp  Gs(k,z+pL/2,z-pL/2)(.  , 


(73-a) 


where 


r r ( , \ ^c2  fzCl-zlL  2 . L(<p)21 

G$  - Fn(K1,<2,PL)  cos  2k  K 8k  J 

x exp{i<  • [z(r]-r£)  + (l-z)(ri-rJ-)  + tv]} 


(73-b) 


and 
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r 


(76-bj 


Gy  * Fn(K1,.c2,pL) 


[z(r|-r^)+(l-z)(rrrj)+T7l| 


Equations  (72)  through  (76)  constitute  the  desired  results  of 
this  section.  In  deriving  these  expressions,  we  have  made  no  assump- 
tions involving  the  size  of  the  outer  scale  of  turbulence,  Lq.  The 
literature  abounds  with  expressions  for  the  log-amplitude  ana  phase 
covariances  and  the  phase-log-amplitude  cross-covariance  for  zero  time 
lag  (t=0)  [58,59]  and/or  axial  receiver  plane  points  (rj=r£=0) 


or  transmitter  plane  points  (r^rj^O)  [22], [61 J, [62].  However t a l 
of  these  derivations  make  use  of  the  assumption  (usually  implicit) 
that  the  outer  scale  is  much  shorter  than  the  range.  We  shall  demon- 
strate that,  for  our  purposes,  no  such  assumption  is  necessary. 


D.  Reduction  of  Expression  for  Spectrum 

In  this  section  the  results  of  the  two  preceding  sections  are 
combined.  Some  simplifying  assumptions  and  observations  are  then 
made  which  tremendously  reduce  the  complexity  of  the  expression  for 
the  scintillation  spectrum. 

Our  procedure  is  to  first  make  the  weak  turbulence  approximation 
(which  we  quantify)  to  simplify  the  expression  for  the  H function  of 
Equation  (29-b).  Then  by  employing  the  expressions  derived  in  the 
previous  section  for  the  required  statistical  moments  and  symbolic- 
ally performing  a number  of  diffraction  integrals  we  arrive  at  the 
final  expression. 

The  H function  of  Equation  (29b)  is  given  by 


(77-a) 


where 


a * 1 [D^iy^sQiO)  +Dw(r3-r4;0;0)] 


(77-b) 


and 
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e s " 1 [DW(rrr4;rrr2;T)  + Dw(r2”r3;ri"r2;T) 

-Dw(rrF3;Fj-F^;T)  - ;T)] 

+ 2[Cx(F1-F3;F]-F^;t)  + Cx(F2-F4;Fj-F^;r) 

+ 1CxS(FTF3;F1-F^t)  “ icxS(F2"F4;Fl'F2;T^  • (77-c) 

Under  the  weak  turbulence  approximation  we  assume  that 

M.  lei  « i (78) 

so  that  we  need  only  keep  the  first  two  terms  of  the  Taylor's  series 
expansions  of  exp(-a)  and  exp(B); 

H % l-(o-e)  - ( 1 — ot ) - e.  (79) 

In  Appendix  C we  delineate  the  conditions  under  which  we  are  justified 
in  making  this  approximation.  For  a von  K^rnrfn  index  of  refraction 
spectrum  these  conditions  are 

a < .546  C2  k2L^Plf3  + p^3)  « 1 

and 


6 < .496  C„  k7/6L11/6  « 1 


where 

Pfj  = 1^,-rjl 


From  the  results  of  the  previous  section  (Equations  (72) 
through  (76))  we  see  that  H (Equation  (29))  is  given  approximately 
by 
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This  expression  can  be  factored  to  yield 


,1/2  f2z  .1  .2(l-z)  , 

H = 2k2L2  jj  dx  / f dz  f dp  + J dz  f dpi  Fn(Kl.K2,pL) 


x e 


'0  J0 
iK-[z(r'rry+Tv] 


’ i (l-z)K*r^-iy  i (l-zjic'^+iy 
e - e 


" -i(l-z)K*r4+iy  -i(l-z)K*r3-iy 
e - e 


(82) 


With  Equation  (82)  plus  the  expression  for  the  spectrum  in  terms 
of  the  H function  (Equation  ( 29-a ) ) we  have 


f00  • 

S (u)  = 2k2L2  dx  e"1u)T 

r —oo 


drjW(rj) 


dr'W(r^) 
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dz 

+ 

dz 

dp 

fi  v 

0 

0 

1/2 

0 

J 

Fn(Kl »K2,pL^ 


x e 


iK*[z(rj-r£)+x7] 


If  we  look  only  at  the  rj,r2  integrals  we  observe  that  they  can 
be  written  in  the  form 


(kl)2  e'iY  |dr,E(r,)  K 


ik  I.?  L(l-z) 

xt-  ^r^-zr-|  • r-j-  v J * 


x dr2E*(r2)  e 


2L  (r2"2r2'rP 


= c/r  , .2L  <rf-2rrri> 


(st)  e'Y  j^EIr,)  e 

x |dF2E*(F2)  e^{r2-2?2-  h * 


t 4.  L0--4  ~ 


I 
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(84) 


Now  recall  the  formula  for  the  Fresnel  diffraction  integral; 

ik  . 7 o—  - 2\ 

->  -ik  ikL  f c/—  \ 2L  (^-2r  -r+r6) 

ri  = p dr$E(rs)  » ^ s s 


E(r) 


2irL 


(85) 


Comparing  Equations  (84)  and  (85)  we  observe  that  the  r-|  and  r2 
Integrations  can  be  performed  symbolically  to  give 


r'y 


-jk 

2L 


_ L(l-z)  -I  - 

1 k 1 


r4rr'2 
>2L  rl 


E*(r,') 


rdJir.2 

e 2L  1 E(rj) 


ik  I-,  A L(l-z)  - 

JT  ri  + “k  K 


E*(ri  * ^T1  «) 


(86) 
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where  it  is  to  be  understood  that  the  fields  are  the  free-space 
receiver  plane  fields.  It  may  easily  be  verified  that  the  r$  and 
r4  integrations  yield,  the  conjugate  of  this  expression  with  the  vari 
able  r-J  replaced  by  r£.  With  these  simplifications  Equation  (83) 
becomes 


SpU)  = 2k2L2 


1 f2(l-z) 

dz  dp 

1/2  JO 
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or 


Sp(«)  = 2k2L2 


2z  fl  f2(l-z) 
dp  + dz  dp  ^ 

0 1/2  Jo 


X ||dtc  Fn(K1,K2,pL)  [2ir6(w-K*v)] 


X |H(z,p,tc,W)  - H*(z,p,-k,W)|2 


(88- 


W (z  ,p  » K »yy)  - 


(88-b) 


In  Appendix  D we  demonstrate  that  due  to  the  scintillation  spec- 
trum's extreme  insensitivity  to  low  spatial  frequency  index  of  refrac- 
tion inhomogeneities,  very  little  error  is  incurred  by  dropping  the 
p dependence  of  the  H function  and  extending  the  limits  on  the  p 
integration  in  Equation  (88-a)  to  infinity  to  yield 


Sp(w)  = 4tt^ k^L  J dz  | J d<  'J>n(K-i,K2»0)6(tl)-K*v) 

^ — oo 

x |H(z,7,w)  - H*(z,-k,w)  | 2 , (89-a) 


where 


H(z,k,w)  - 


drW(r)elK’r  E^r 


(89-b) 


and  we  have  made  use  of  the  relationship  between  the  two  and  three 
dimensional  index  of  refraction  spectra, 


dp  Fn(Kl  »k2’pL)  = f *n^iel*ie2*0^ 

J -oo  L 


Finally  without  loss  of  generality  we  can  assume  the  wind  to  be 
blowing  in  the  x direction, 

7 = x v , (91) 

so  that  performing  the  q integration  in  Equation  (89)  yields 


Sp(o>)  = (4*2k2L/v) 


,1 

dz 

0 


„d*  *„  (a-  ,k,o) 


x |H(z,k',w)  - H*(z,-k'  ,W) 1 2 


(92) 


where 


and  we  have  dropped  the  subscript  on  the  remaining  spatial  frequency 
variable  of  integration. 

Equation  (92)  is  the  expressed  objective  of  this  chapter.  The 
utility  of  this  formula  lies  in  the  fact  that  it  expresses  the  temporal 
scintillation  spectrum  for  an  arbitrary  field  and  that  the  expression 
for  the  source  field  enters  the  formula  only  in  terms  of  the  free- 
space  receiver  plane  fields.  This  latter  fact  plays  an  important  ro^ 
in  the  chapters  to  follow. 

In  addition  to  the  mathematical  simplicity  of  Equation  (92), 
its  various  components  are  easily  interpretable  in  terms  of  a physical 
model  of  the  turbulent  atmosphere.  The  H function  in  particular  dis- 
plays some  very  interesting  behavior  which  offers  insight  into  the 
propagation  problem.  Equation  (92)  is  explored  and  discussed  more 
fully  in  the  next  section. 


E.  Physical  Model 

Now  we  wish  to  interpret  our  expression  for  the  temporal  scin- 
tillation spectrum  in  the  context  of  a physical  model  of  the  atmo- 
spheric turbulence.  In  this  model  the  index  of  refraction  inhomo- 
geneities are  envisioned  as  a series  of  phase  gratings  [63]  of 
various  orientations  and  spatial  frequencies.  The  components  of 
Equation  (92)  will  then  be  shown  'to  describe  the  behavior  of  the 
fields  in  terms  of  the  properties  of  these  phase  gratings. 

We  wish  to  interpret  Equation  (92)  of  the  previous  section: 


x |H(z,k' ,w)  - H*(z,-k' ,w) I2 


(93) 
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where 


H(z,k',w)  = 

-i (L-z)k1^  _ _ 

e 2k  |drW(r)eiK  ’r  E^r  - E*(r)  , 

«•  ■ (7  •")  • <94> 

and  we  have  changed  the  z variable  of  integration  to  the  product  L*z. 

As  preliminary  background  to  this  discussion  we  introduce  the 
concept  of  a phase  grating  as  sketched  in  Figure  7.  Illustrated  in 
this  figure  is  an  arbitrary  wavefront  impinging  on  a phase  grating. 
The  phase  grating  is  merely  a diffraction  grating  whose  transmission 
function  is  periodic  in  phase  rather  than  amplitude.  To  the  right  of 
this  phase  grating  (which  can  be  thought  of  as  a periodic  variation 
in  index  of  refraction  along  a plane  perpendicular  to  the  direction 
of  propagation)  we  have  depicted  only  the  zeroth  order  and  plus  and 
minus  one  order  diffracted  wavefronts.  It  is  easily  shown  [64]  that 
the  angles  of  the  axes  of  the  two  diffracted  wavefronts  with  respect 
to  the  axis  of  the  undiffracted  wavefront  are  given  by 

a = x/i  , (95) 

where  A is  the  wavelength  of  the  field  and  1 is  the  period  of  the 
grating. 

Now  consider  the  situation  shown  in  Figure  8.  Here  we  have  a 
phase  grating  of  period 

2ir/ k (96) 

and  orientation 

k/M 

located  in  the  plane  z. 


An  axial  ray  striking  this  phase  grating  from  the  left  will  produce 
zeroth  and  plus  and  minus  one  order  diffracted  rays  to  the  right 
From  the  previous  discussion  the  diffraction  angles  are  seen  to  be 


and  in  the  plane  z=L  the  displacements  from  the  propagation  axis  are 


(98) 


But  this  is  merely  the  displacement  of  the  field  in  the  expression 
for  the  H function. 


The  foregoing  argument  was  an  heuristic  derivation  of  the  field 
translations  observed  in  the  H function.  In  Appendix  D it  is  shown 
much  more  rigorously  that  the  total  diffracted  field  of  a weak 
sinusoidal  phase  grating  is  given  by 


Et  = eik(L-z)  | c E 
n=-l  n n 


(99-a) 


where 


En  = e 


_ .micT 


-in2(L-z)  2 

2k  K 


E(F-=^7)  , 


(99-b) 


'n  * 0 f) ' " ' 


(99-c) 


and  c is  a small  number  proportional  to  the  "strength"  (maximum  phase 
change)  of  the  phase  grating  In  the  above  equations,  n represents 
the  diffracted  order  and  E(r)  is  the  incident  field. 


The  total  irradiance  is  given  by 

,2.  1 1 


’T  ■ IM  - I I VS  EnEi!  • 
n=- 1 m — I 

and  the  total  receiver  power  is 
PT  = | dF  W(f)  IT 


O00) 


1 1 


J.,  j-i  c-c;  I M(F)  E-E" 


(101) 


L 
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We  retain  only  the  terms  in  e°  and  e‘  (because  e is  small)  so  that 
Equation  (1C1)  gives 


Py  = dr  W(r)  E(r)  2 
* 


+ C-1C0  e 


-HL-Z)  k2 

2k  K 


+ c-ico  6 


lilyil  k2 

2k 


+ CoC*  e 


IjLrJd  2 

2k  K 


+ C0C1  e 


KLzzl  2 
2k 


r 


dr  W(r)  e lK  r 


| dF  W(r)  el  K r 
I dr  W(F)  e"ilc  r 
| dr  W(r)  ei,c*r 


E(r  + 1Lrl  *)  E*{r) 

E*^r  + K^j  E(r) 

E*(F  - 7^  E(r) 

- ^L~z^  7j  E*(r)  . 

002) 


Each  of  these  integrals  (with  the  exception  of  the  first)  is  merely 
the  H function  or  its  complex  conjugate.  The  H function  is  therefore 
seen  to  represent  the  interference  between  the  zeroth  order  field  and 
the  plus  or  minus  one  order  field. 

With  the  utility  of  the  phase  grating  model  of  the  turbulence 
field  thus  demonstrated,  we  shall  now  demonstrate  that  the  general 
scintillation  spectrum  formula  (Equation  (93))  can  be  derived  using 
this  model  plus  a few  previously  adopted  definitions. 

Essentially  what  we  have  implied  in  the  preceding  development  is 
that  the  index  of  refraction  field  can  be  expanded  in  an  orthogonal 
set  of  (moving)  phase  screens,  the  transmission  functions  of  which 
are  given  by  exp(ie  cosK*r).  The  problem  however,  is  that  the  set  of 
cosine  functions  in  itself  is  not  complete  [65].  Since  the  functions 
sine  plus  cosine  do  compose  a complete  orthogonal  set  [66],  we  now 
wish  to  apply  the  preceding  technique  to  the  situation  in  which  the 
"phase"  grating  transmission  function  is  given  by 


where  e is  a small  (generally  complex)  number.  From  Appendix  E,  the 
resulting  receiver  plane  field  is  seen  to  be 


. -in2(L-z)tc2 

ET(z,7,t)  * e,k(L'z)  l cne,n‘-re 

n=0 

E(F‘n^ri'7)  * 004-a) 

where 


C0  = 1 , (104-b) 


Cj  = ie(z,K,t)  , (104-c) 


and  we  have  indicated  explicitly  the  dependence  of  the  field  upon  the 
path  variable,  z,  the  spatial  frequency,  «,  and  the  time,  t.  In  the 
strength  function,  £,  the  z dependence  denotes  the  location  of  the 
phase  grating,  the  < dependence  indicates  the  strength  to  be  a function 
of  spatial  frequency,  and  the  time  dependence  reflects  the  fact  that 
the  phase  grating  is  being  blown  by  the^wind  across  the  propagation 
path. 


Following  arguments  and  definitions  employed  earlier  in  this 
section.  Equations  (94),  (102),  and  (104)  give  for  the  receiver  plane 
power 


PT(z,K,t)  = Pnc  - ie*(z,K,t)H*(z,K,w)  + ie(z,K,t)H(z,K,w)  , 

(105-a) 


where 


PDC  5 lcol2  } dr  W(r)|E(r)|2 


(105-b) 


As  in  Section  B of  this  chapter,  we  define  the  scintillation  spectrum 
as  the  Fourier  transform  of  the  power  covariance 


(106) 
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where  the  power  covariance  is  given  by 


Cp(x)  = f^dz-j  J^dz2  jdw  JdK*  <[PT(z1,<,t)-PDC] 


x [Pyv z2 » K ' »t+T)-PDC^*> 


007) 


In  this  equation  we  have  "summed"  the  power  (density)  of  Equation 
(105)  over  all  spatial  frequencies  and  over  all  phase  grating 
locations  between  zero  and  L. 


With  Equations  (105)  and  (107)  the  power  covariance  becomes 

- tdii  j*t  F K 

X [<e (z-|  ,K,t)e*(z2.K'  • t+T ) >ff ( Z 1 ,k,w)H*(Z2,k  1 ,w) 


- <e(z1,K,t)  (z2,k' ,t+T)>H(z-|  ,k,w)H(z2»k' ,w) 


- <e*(z1  ,ic,t)e*(z2,ic'  ,t+r)>H*(Zi  ,k,w)H*(z2,k' ,w) 


+ <e*(z1 ,K,t)e(z2,K' ,t+x)>H*(z1 ,k,w)H(z2,k' ,w)] 


(108) 


We  now  claim  that  the  correlation  function  of  the  phase  grating 
strength  is  given  by  [67]  (see  Equation  (58)) 


<e(z1,K,t)e*(z2,K' ,t+t)>  = k2  eiTIC'*V  6(|k-k'|) 


x Fn^Kl  ,K2 » I Z1  “z2"l ) * 


(109) 


where  F0  is  the  two  dimensional  spatial  spectrum  of  the  index  of 
refraction  field.  In  this  equation  the  exponential  phase  term 
reflects  the  frozen  turbulence  approximation  and  the  statistical 
orthogonality  (with  respect  to  spatial  frequency)  is  indicated  by 
the  delta  function.  This  orthogonality  states  that,  on  the  average, 
only  phase  grating  pairs  with  equal  spatial  frequencies  and  identical 
orientations  contribute  to  the  fluctuations  in  the  received  power. 
Finally,  the  square  of  the  wavenumber,  k,  yields  the  correct  units. 
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We  also  claim  that,  for  the  term  in  Equation  (108)  in  which  neither 
of  the  strength  functions  is  conjugated,  the  correlation  function  is 
given  by  [68] 


<e(z-|  ,K,t)e(z2,K'  ,t+r)>  = k2  e"iTK  *V  6 ( | k+k  ' | ) 

x Fn(Kl  ,K2’  I Z1  -z2 1 ) • (no) 

Combining  Equations  ( 1 08 ) » (109),  and  (110)  and  performing  the  7' 
integrations  yields  for  the  power  covariance 


Cp(t)  = k2  | dz-,  J dz2  |d<  Fn(K1,K2,|z1-z2|) 


x < e 


t 


0 

.1  TK*  V 


- e 


[H( z-j  ,k,w)H*(z2,k,w)-H(z-|  ,k,w)H(z2,k,w)] 

"1tic*v  [H*(z-j  ,7,w)H*(z2,-7,w)-H*(z-]  ,7,w)H(z2,7,w)]|  . 


(Ill) 


Assuming  isotropic  turbulence  and  making  use  of  the  change  of 
variable 

K *■  -K 

gives 

cp(x)  = k2  | dzi  JQdz2  |d7  Fn(<1  ,tc2.  |zi-z2|  )elT<’v 
x [H(z1 ,7,w)H*(z2,7,w)-H(z1,7,w)H(z2,-7,w) 


- H*(z-j,-K,w)H*(z2,-K,w)+f/*(z-j,-ic,w)H(z2,-K,w)]  . (112) 


Switching  to  sum  and  difference  coordinates  on  the  range  variables 

z = (z-|+z2)/2  , p = zrz2 

(which  correspond  respectively  to  the  mean  position  and  separation 
of  the  phase  grating  pairs)  and  using  the  argument  that  the  index 
of  refraction  spectrum,  Fn,  is  very  small  for 


(see  the  discussion  in  Appendix  D)  we  obtain  finally 

Cp(f)  = 2nk2  J^dz  Jd<  eiTK  V (^i 

x |H(z,k,w)-H*(z,-iT,w)|2  . (113) 


The  scintillation  spectrum  is  then  given  by 

L OO  X 

Sp(u))  = (4ir2k2/v)  j^dz  | dK  *n  ,K,oj 


x |H(z,k',w)-H*(z,-(c',w)|2  , (114) 


where 

k1  = (w/v.k)  , 

and  we  have  assumed  the  wind  to  be  blowing  in  the  x-direction. 
Comparison  of  Equations  (93)  and  (114)  shows  them  to  be  identical. 

It  is  indeed  satisfying  that  with  a few  reasonable  assumptions  we 
have  duplicated  the  results  of  the  much  more  rigorous  development 
of  the  preceding  sections. 

The  purpose  of  this  section  has  been  to  give  a more  physical 
picture  of  the  mechanisms  giving  rise  to  temporal  scintillation. 
Salient  features  of  the  development  contained  herein  were  the  expan- 
sion of  the  index  of  refraction  field  in  a complete  orthogonal  set 
of  phase  screens  and  the  symbolic  performance  of  the  diffraction 
integral  to  obtain  the  fields  diffracted  by  the  phase  screens.  By 
using  previously  adopted  definitions  for  i rradiance,  power,  and  the 
temporal  scintillation  spectrum,  and  making  some  physically  justi- 
fiable assumptions,  the  results  of  the  preceding  sections  were 
duplicated.  The  H function  in  particular  was  seen  to  have  an  inter- 
esting interpretation  in  the  context  of  a phase  grating  model  of 
the  atmosphere.  It  represents  the  power  associated  with  the  inter- 
action between  the  zeroth  and  plus  or  minus  one  order  diffracted 
fields  of  the  phase  grating. 

This  concludes  the  discussion  of  the  formula  for  the  temporal 
scintillation  spectrum. 


51 


r 1 


F.  Summary 


In  this  chapter  we  have  performed  the  bulk  of  the  mathematical 
development  which  is  to  be  applied  in  subsequent  chapters.  The 
objective  was  to  derive  a general  expression  for  the  temporal  scin- 
tillation spectrum  of  an  arbitrary  source  field  and  an  extended 
receiver.  This  was  accomplished  in  several  stages.  Use  was  made  of 
the  Extended  Huygens-Fresnel  integral  to  derive  a formula  for  the 
spectrum  in  terms  of  the  second  order  statistical  moments  of  spheri- 
cal waves.  Expressions  for  these  statistical  moments  (within  the 
weak  turbulence  regime)  were  derived  by  the  method  of  smooth  pertur- 
bations and  applied  to  this  formula.  Various  simplifications  and 
manipulations  finally  resulted  in  an  extremely  compact  yet  versatile 
expression  which  was  then  interpreted  using  a phase  grating  model  of 
the  turbulent  atmosphere. 

Section  B of  this  chapter  was  devoted  to  developing  a general 
expression  for  the  temporal  scintillation  spectrum  of  an  arbitrary 
field  source.  An  extension  of  the  familiar  Huygens-Fresnel  diffrac- 
tion integral  was  employed  to  express  the  receiver  plane  field  of  an 
arbitrary  source  in  terms  of  the  complex  phase  perturbations  to  a 
spherical  wave.  The  power  received  at  the  detector  wr  defined  as  the 
surface  integral  of  the  irradiance  (square  modulus)  o^  ..iis  field,  and 
the  temporal  scintillation  spectrum  as  the  Fourier  transform  of  the 
covariance  of  this  power.  The  resulting  formula  required  knowledge  of 
the  covariance  of  irradiance  of  spherical  waves.  Under  the  assumption 
that  the  spherical  wave  complex  phase  perturbations  were  normally  dis- 
tributed, this  fourth  moment  of  the  field  was  expressed  as  a series  of 
second  order  statistical  moments.  One  interesting  aspect  of  the 
general  result  of  this  section  was  that  the  expression  for  the  spectrum 
retained  terms  proportional  to  the  cross-covariance  of  phase  and  log- 
amp  1 itude. 

Expressions  for  the  required  log-amplitude  covariance,  phase-log- 
amplitude  cross-covariances,  and  wave  structure  functions  were  derived 
in  Section  C via  the  method  of  smooth  perturbations.  These  derivations 
relied  upon  the  assumption  of  local  homogeneity  of  the  index  of 
refraction  statistics,  and  the  paraxial  ray  approximation.  Temporal 
behavior  of  the  statistical  moments  was  deduced  by  making  use  of  the 
frozen  turbulence  hypothesis.  The  results,  of  this  section  showed 
the  desired  statistical  moments  to  be  simultaneously  homogeneous  in 
temporal  separations  and  spatial  separations  within  both  the  trans- 
mitter and  receiver  planes. 

The  results  of  Sections  B and  C were  combined  in  Section  D.  By 
using  the  weak  turbulence  approximation  and  making  some  simple  alge- 
braic manipulations,  the  resulting  formula  was  expressed  in  terms  of 
a series  of  diffraction  integrals  which  were  performed  symbolically 
to  yield  the  final  expression  for  the  spectrum.  This  formula  for  the 
temporal  scintillation  spectrum  of  an  arbitrary  source,  was  in  terms 
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of  integrations  over  the  receiver  plane  fields,  a single  spatial  fre- 
quency variable,  and  a propagation  path  variable. 

In  Section  E a model  of  the  atmospheric  turbulence  was  developed 
in  terms  of  weak  phase  gratings.  Using  diffraction  theory  and  some 
physically  justifiable  assumptions  the  results  of  the  much  more 
rigorous  development  of  the  preceding  sections  were  duplicated,  and 
the  component  parts  of  the  formula  for  the  scintillation  spectrum 
were  given  physical  interpretations. 

The  following  chapters  are  concerned  with  the  application  of  the 
results  of  this  chapter  to  problems  of  contemporary  interest. 
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CHAPTER  III 

ANALYSIS  OF  A CLEAN  GAUSSIAN  BEAM 


A.  Introduction 

This  chapter  is  concerned  with  application  of  the  results  of 
Chapter  II  to  calculation  of  the  temporal  spectrum  of  a gaussian 
beam  with  an  off-axis  point  detector.  The  choice  of  this  particular 
source  description  is  appropriate  because  the  output  of  most  lasers 
is  closely  approximated  by  a gaussian  beam.  Further,  the  assumption 
of  an  off-axis  detector  is  reasonable  if  for  instance  the  laser  beam 
is  being  steered  by  a servo  system.  Nominally  a tracker  would  steer 
the  laser  beam  so  that  it  is  centered  upon  the  detector.  Since  how- 
ever this  may  not  always  be  the  case  [69],  the  development  herein 
is  of  obvious  importance. 

Within  this  chapter  an  evaluation  of  the  temporal  scintillation 
spectrum  of  a gaussian  beam  in  the  presence  of  weak  turbulence  and 
with  an  off-axis  detector  is  presented.  Analytic  expressions  are 
developed  for  the  asymptotes  of  the  spectrum  with  axial  as  well  as 
off-axis  point  detectors.  The  orientation  of  the  off-axis  detector 
with  respect  to  wind  direction  is  shown  to  have  important  effects 
within  the  low  frequency  portion  of  the  spectrum.  For  a typical 
laser  beam,  computer  analyses  of  the  spectrum  are  provided.  Plots 
are  presented  of  the  differential  path  contributions  as  well  as  the 
spectra  for  various  detector  locations.  Finally  we  discuss  the  ex- 
tension of  the  theory  presented  herein  to  the  treatment  of  finite 
receiver  apertures  and  a finite  inner  scale. 

As  a result  of  the  analysis  it  is  shown  that  for  an  off-axis 
detector,  the  low-frequency  behavior  of  the  gaussian  beam  spectrum 
departs  markedly  from  that  of  plane  and  spherical  waves.  In  ad- 
dition, the  differential  path  contribution  for  an  off-axis  detector 
is  peaked  toward  the  transmitter  (as  it  is  for  a plane  wave)  even 
though  the  receiver  plane  is  well  within  the  far  field  of  the  laser 
beam. 

In  section  B of  this  chapter  use  is  made  of  the  general  formula, 
derived  in  Chapter  II,  in  developing  an  expression  for  the  temporal 
scintillation  spectrum  of  a gaussian  beam  with  an  arbitrarily  placed 
point  detector. 

An  asymptotic  evaluation  of  the  gaussian  beam  scintillation 
spectrum  is  given  in  section  C.  Analytic  expressions  for  the 
asymptotic  high  and  low  frequency  behavior  are  provided. 


Section  D contains  numerically  calculated  spectra  for  several 
typical  beam-detector  orientations.  In  this  section  we  also 
present  plots  of  the  computer-calculated  differential  path  con- 
tributions for  some  typical  cases. 

Two  additional  topics,  the  effect  of  finite  receiver  apertures 
and  non-zero  inner  scales,  are  discussed  in  section  E. 

Section  F summarizes  the  chapter. 

B.  Development  of  Expression  for  Gaussian  Beam  Spectrum 

We  now  wish  to  apply  the  results  of  Chapter  II  to  derive  an 
expression  for  the  scintillation  spectrum  of  a gaussian  beam. 

Throughout  this  section  and  the  remainder  of  this  work  we 
shall  assume  that  a laser  beam  field  is  typified  by  a spherical 
phase  front  and  a gaussian  field  amplitude  profile  [39] 

, ik(x2+y2)  _ (x2+y2)T 
E(x,y,z)  = Eq(z)  el  2R(z)  w2(z)  J , 


where  R(z)  and  w(z)  are  real  functions  describing  respectively  the 
radius  of  curvature  and  the  e_l  amplitude  (or  spotsize)  of  the  beam. 
To  simplify  notation  however  we  will  write  the  field  in  terms  of  a 
complex  radius  of  curvature 


-k(x2+y2j 


E(x,y,z)  = EQ(z)  e 


Specifically  we  assume  that  the  gaussian  beam  field  within  the 
receiver  plane  is  0 

-kN/2L  rc 

E(r)  - e1kUltVL)  Jiia  J+iNinyu 


l+iN(l+Z  /L] 


(115a) 


where  the  beam  waist  is  located  at  z=-z  and  N is  a Fresnel  number 
defined  as 


N 


(115b) 
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where  w0  is  the  beam  waist  spotsize.  For  this  field  the  receiver 
plane  radius  of  curvature  and  beam  spotsize  are  respectively 


and 


R(L)  « L 


1+N2(1+z0/L) 

N2(1+z0/L) 


W(L)  = wQ  Jl+N2(l+z0/L)2. 

For  a point  detector  located  at  rQ; 

w(r)  = 6(|r-r0|) 

the  H function  of  Chapter  II  becomes 


H(z,k',w)  = 


(I-z)Lk'2 

2k 


k 


E*(r0)  . 


Use  of  Equation  (115a)  for  the  receiver  plane  field  gives 


where  we  have  defined 


with  some  straightforward  algebra  we  see  that 


I H{z , K ' *W ) - ^(z.-ic1  ,w)  r = 

-(2kN/L)r2  -N(l-z)2L/k  ,2 

2 1+N2(Ue)2  l+N2(I+e)2 


kN/2L 


l+N2(l+e)2 


x “S  2 cosh 


2N(l-z) 

Ll+N2(l+e)2  '° 
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(119) 
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'i 

l+N2(l+e) (z+e)  r 

k K 

- 1+N2(1+e)2  J 

By  employing  Equation  (92)  of  Chapter  II  and  the  double  angle 
formulae  for  the  trigonometric  and  hyperbolic  functions  we  then 
obtain  the  expression  for  the  gaussian  beam  scintillation  spectrum 


S (u.)  = (16ir2k2L/v)  pi 

P Ll+N2(l+e)2J 


-(2kN/L)r^ 

l+N2(l+e)2 


-N(l-z)^L/k  .2 


' dz  f dK.nft  ,.,o) 

n ' -oo  Vv  ' 


(120) 


57 


For  the  special  case  of  an  on-axis  point  detector  (ro=0),  Equation 
(120)  reduces  to 


Sp(«o)  = (327r2k2L/v) 
-N( 1 -z)2  .2 


kN/2L 


Lhn2(i+c)2J  Jo 


dz 


dK»„(f«,0) 


:HN2(H£)2  s1n2)li^lkK.2 


l+N2(l+e)(z+e) 

i , n2/ 1 _ \2 


1+fTO+e) 


(121) 


Equation  (121)  was  derived  and  evaluated  for  a variety  of  Fresnel 
numbers  and  beam  waist  locations  by  Ishimaru  [39].  He  however 
employed  the  technique  of  Tatarski  [70]  in  obtaining  this  expression 
for  the  spectrum.  It  is  pleasing  that,  under  certain  circumstances, 
these  two  techniques  (Tatarski 's  and  the  Extended  Huygen' s-Fresnel ) 
should  produce  exactly  the  same  results. 

In  the  limit  as  N+0  (plane  wave)  and  N->°°  (spherical  wave)  Equation 
(120)  yields  respectively 


sp(w) 


1 

dz 

0 


dK<(>n  ,K,oj  sin2  ^ k,2J  (122a) 


and 


Sp(“) 


dz 


dK*n  (n  ’K’°)  slV 
o'7 


(l-z)L  z±£  k.2*1 

2k  1+e  J 


(122b) 


Equations  (122a)  and  (122b)  which  are  expressions  for  the  plane  and 
spherical  wave  scintillation  spectra,  are  treated  extensively  by 
Tatarski  [36]  and  Clifford  [35]. 

Throughout  the  remainder  of  this  work  use  will  be  made  of  the 
von  Karman  form  of  the  three-dimensional  index  of  refraction  spatial 
spectrum  [31], 


<f,n(Ki  »k2,k3J  = -°33C 


2,  2,  2 . 

k1+k2+k3  + 


(^J 


-11/6 


(123) 


whereas  most  previous  workers  [35],  [39],  [71]  assumed  the  Kolmogorov 
(infinite  outer  scale)  model.  In  addition.  Equation  (120)  can  be 
put  in  a form  more  convenient  for  analysis  by  non-dimensionalizing  the 
•c  integral  through  the  change  of  variable 
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K'  • (124) 


The  temporal  spectrum  is  then  given  by 


where  all  of  the  primed  variables  denote  those  which  have  been 
normalized  to  the  plane  wave  Fresnel  zone: 


o 

-J 

II 

- o 

-J 

(125b) 

K1  = /aLk  = 2tt/(  A//TT) 

(125c) 

r'  = rJ^L 

0 0 

(125d) 

and 

v/JIl  s o)q  . 

(125e) 

Since  the  predominant  diffraction  pattern  size  at  the 
dimension  on  the  order  of  the  Fresnel  zone,  711  [62], 

receiver  is  of 
the  normalization 
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J 


frequency,  u0,  is  seen  to  be  the  temporal  frequency  observed  when  this 
scale  size  is  blown  through  the  beam  path  at  velocity  v. 


We  now  write  Equation  (125a)  in  the  form 
SD(w)  « .033C^47T2k2L(v’AL)8/3v-1  k 

P L l+NZ(l+e)2 

-4*Nr'2 
o 

x el+N2(l+e)2  CS1  (u>)+S2(co)], 


where  the  functions  S^w)  and  S2 (to ) are  given  by 


5^  (w)  = 4 


-N(l-z) 


: -L-i-’ tHT/i 


2n[l+N2(l+£)2] 


.2 J N(l-Z 


(126a) 


(126b) 


,(  x‘  SL.  +y 1 k‘ 

1+N2(1+e)2  ' 0W0  0 


S2(u>)  =4  dz  d< 1 1 < 1 
o -°°  t 


-N(l-z)‘ 


* OP 
■'Hkfl 


(126c) 


2tt[1+NZ(  l+e)Z] 


l+N2(l+e)2  JL 


%)}}■ 


The  purpose  of  this  division  of  the  spectrum  into  two  components  is 
to  isolate  the  effects  due  to  the  detector  being  off-axis  (given 
by  the  function  S](w))  from  the  behavior  of  the  spectrum  with  the 
on-axis  detector  (given  by  the  function  S2(w)).  If  the  spectrum 
of  Equation  (126)  is  normalized  to  the  spherical  wave  log-amplitude 
variance  [21] 


2 = .Q99tt7/2  21/3 
S 64r(4/3)(-cos  tt ) 


TT— CV'V1'6 
II  \ n 


(127) 
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we  obtain  finally 


u,0s  (u>) 

“5Ji 7 

r(0.0,L)of 


.,2 


-4uNr0 

..2/, . \Z 


= Ce'+N(,+e)  [S,(«)«2(.)] 


(128a) 


where  I is  the  on-axis  receiver  plane  irradiance 
\/  kN/2L 


1(0, 0,L)  = 


l+iN(l+e) 


(128b) 


and 


c = 156_  ^ r(4/3)(-cos  w)  = 48.527.  (128c) 

3it 

This  particular  normalization  is  chosen  so  that  the  area  under  the 
spectrum  is  normalized  to  the  spherical  wave  log-amplitude  variance; 

r°°  u)  S (w) 

lim  d -8-E T = 4 x 2™,  (129) 

J-»  r(0,0,L)o‘ 

e-0 

I -x30 

0 

where  the  additional  factor  of  four  arises  from  the  difference  between 
the  variance  of  irradiance  and  the  log-amplitude  variance. 


under  the  weak  turbulence  approximation  [62].  Equation  (129)  can  be 
verified  by  taking  the  indicated  limits  of  expression  (129),  changing 
the  k.'  ,w  variables  to  polar  coordinates,  performing  the  angular  in- 
tegration, and  employing  integration  formulae  presented  by  Tatarski 
[72]. 

Equation  (128a)  is  in  a form  convenient  for  asymptotic  and  numer- 
ical analysis.  The  asymptotic  evaluation  for  each  of  the  two  com- 
ponent spectra  (Si  and  S2)  of  this  expression  will  be  presented  in 
the  next  section. 


61 


C.  Asymptotic  Evaluation 

In  this  section  we  develop  analytic  expressions  which  describe 
the  asymptotic  high  and  low  frequency  behavior  of  the  scintillation 
spectrum  of  a gaussian  beam  with  an  off-axis  detector.  These  formu- 
lae are  derived  assuming  an  arbitrary  beam  Fresnel  number,  N,  and 
an  arbitrary  detector  location,  r0.  Because  the  spectra  Si(w)  and 
S2(w)  exhibit  diverse  behaviors,  they  are  evaluated  separately. 

Asymptotic  Behavior  of  S^(oj) 

The  Si  spectrum  of  Equation  (128a)  arises  simply  because  the  de- 
tector is  off-axis  to  the  gaussian  beam.  It  is  identically  zero  if 
the  detector  is  axially  located.  Typical  results  presented  later  in 
this  section  will  show  that  Si  makes  a significant  contribution  to  the 
total  spectrum  only  within  the  low  frequency  range,  i.e.,  for  fre- 
quencies such  that 


— < 2tt. 

(i)  'V* 


(130) 


We  now  wish  to  make  two  approximations  which  will  greatly  simplify 
tne  evaluation  of  the  spectrum  Si.  First  we  approximate  the  exponential 
function  and  then  the  hyperbolic  sine  function  in  Equation  (126b).  The 
argument  of  the  exponential  function  is  bounded  in  its  z and  N variables 
in  the  following  manner: 


N 1-Zr  r i 2 . / , ,2-,  1 r ,2.,  , x2n 

J [k  +(u/w  ) ] < 37  [k  + (uj/u)  ) ]. 

27T[l+N2(l+e)Z]  0 * 0 


At  the  point 


<•  = 7(^/^0)2  + 0.077/L^)2, 


beyond  which  the  function 


[k|2+(w/w  )2+(l . 077/L ' )2]— 1 1/6 


decreases  rapidly,  we  have 


031) 


(132) 


(133) 


h [k,2+(«/»0)2]  = 37  2(®/«0)z+  ( 


034) 


If  we  restrict  our  attention  to  temporal  frequencies  such  that 


Lq  > 1.077 


201 

4tt 


AL  <x.  4j  XL 


we  can  use  the  approximation 


0 37) 


-N(1-z)2[k,2+(u/uo)2] 

e27r[l+N2(l+e)2] 


(138) 


We  now  must  approximate  the  hyperbolic  sine  function  in  Equation 
(126b).  Its  argument  is  bounded  as  follows: 


NQ-z) 


l+N2(l+e)2 


(x'  — + y' k'  ) < 

0 wo  0 ^ 2AL 


(x„  — + yrt  k'). 

0 coQ  •'O 


(139) 


At  the  point 

■c’  -P.f  + (1.077/L^)1 
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we  find 


E(xn‘V)s  fa  <~/»o>2+<’-°7'/l;>2> 


(140) 


For  temporal  frequencies  such  that 

- 4 '0  ft?2) 

“o'  ' Lo  / 

the  argument  of  the  sinh  function  is 

N_P-z) /„.  < o_A  5 (1 .077WTL  , v 

i+N2(i+e)z  ' ° «.  V > i nr  "“LT  (V*.» 


< 5(1.077) 

1 Lo  °' 


(141) 


By  requiring 
5(1 . 077)r 


Lo  « 1 


or 


ro  <<5(l .077) ’ 

we  can  employ  the  approximation 
sinh2  X % X2  . 

With  the  foregoing  approximations  the  expression  for  Si  in 
Equation  (126b)  becomes 


(142) 


(143) 


slU)  * j 


1+N2( 1+e )2 


-11/6 


y;*1)2. 


O (0Q  '0 


(144) 


By  expanding  the  quadratic  term  in  the  integrand  of  Equation  (144)  it 
is  seen  that  the  cross  term  is  an  odd  function  of  and  the  integral 
goes  to  zero.  Therefore  Equation  (144)  is 


v«)  = — jr -7  1 

1 JLl+r(l+e)  J 


V 

+ (oj/oJq)^  h 

' 0 - 

j2j-n/6 

(145) 

CO 

+ (y^)2 j <' 

—00 

2|k'2+(<o/o)0)2  +( 

r1.077f 

rcr) 

]-n/6. 

The  change  of  variable 

■ ■ -iii-j  (m2 f\ 


allows  the  expression  of  these  integrals  in  terms  of  Beta  functions  [73]; 


S-j  (w) 


4 T N 

3 [l+N2(l+eV 


{*?  ala'  • raT] 

-4/3 

5 *4$  • «?]""• 

’ , 047) 

where  we  have  used  the  well  known  identity  relating  the  gamma  (r)  and 
beta  functions  [73]. 
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For  the  detector  and  beam  axis  along  a line  col  linear  with  the 
wind,  rg  = (x^,0),  (recall  that  we  assumed  the  wind  to  be  blowing  in 
the  x direction)  then  Equation  (147)  predicts 


S 


1 


a 


u»/w0  < 1 .077/Lg 
w/w0  > 1 .077/1-^ 


(148a) 


However  for  the  detector  and  beam  axis  along  a line  orthogonal  to  the 
wind  direction,  r^  = (0,y^), 


constant 

(u)/w0)"2/3 


; (o/u0  < 1 - 077/ Lq 
; u/u  > 1.077/L^ 


(148b) 


These  asymptotic  results  are  summarized  in  Figure  9a  in  which  we  have 

defined  a as 
r 


0r  l+N2(l+e)2 
and  the  constant  C-j  as 

c _ ^r(4/3) 

1 r(ll/6). 


(149a) 


(149b) 


Figure  9b  shows  the  resulting  numerical  values  of  the  asymptotes  when 
N = 2 x 103 

e = 0 (150) 

L = 10  meters 
o 

JIT  = 2. 51 6x1 0"2  meters 
and  the  detector  is  located  at 
F0  = (.63,0) 
and 

rQ  = (0,-63) . 


(151a) 

(151b) 
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This  can  be  verified  by  differentiating  Equation  (147)  with  respect 
to  w/iOg. 

We  now  return  to  our  phase  grating  model  of  the  atmospheric 
turbulence  and  sketch,  as  in  Figure  10a,  the  receiver  plane  positions 
of  the  zeroth  and  plus  and  minus  one  order  diffracted  beams  as 
functions  of  the  spatial  frequencies  u>/v  and  k.  (Recall  that  the 
spatial  frequency  w/v  is  associated  with  the  wind  direction.)  This 
figure  merely  illustrates  the  fact  that,  for  a finite  beam,  the 
detector  position  with  respect  to  the  wind  direction  and  the  zeroth 
order  beam  axis  will  have  an  effect  upon  the  scintillation  spectrum. 
Of  course  for  fields  of  infinite  extent,  e.g.,  plane  and  spherical 
waves,  the  detector  "position"  loses  distinction.  Indeed,  it  can 
easily  be  seen  that  the  S-j  (oj)  spectrum  goes  to  zero  as  the  beam 
Fresnel  number  approaches  zero  or  infinity.  These  limits  correspond 
respectively  to  plane  and  spherical  waves. 
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LOCUS  OF  *1  ORDER 
OIFFRACTED  BEAM 
CENTER  AS  *-*>00 


For  the  detector  off^axis  to  the  zeroth  order  beam  the  0-  1 
beam  interaction  and  the  1-0  beam  interaction  contribute  un- 
equally to  the  scintillation  spectrum.  When  the  detector  is  in 
the  Is*  or  4th  quadrants  the  contribution  to  the  spectrum  is 
weighted  toward  the  0-+l  beam  interaction.  Likewise,  when  the 
detector  is  in  the  2nd  0r  3rd  quadrants  the  "1-0  beam  inter- 
action predominates. 

In  Figure  10b  we  have  illustrated  the  zeroth  and  +1  order  beams 
for  a phase  screen  located  at  z=0.  Also  shown  in  this  figure  is  the 
beam  deflection  region  for  which  the  index  of  refraction  spectrum, 
4>n,  is  constant 

I 2 , , , \2  L L 1.077 
+ («/v)  k £ k — r-  . 


With  this  figure  we  wish  to  illustrate  the  behavior  of  the  spectrum 
for  the  case  in  which  the  detector  is  at  either  zero  or  90  degrees 
with  respect  to  the  x-axis. 

For  the  zero  degree  case  (ro=(xo,0))  the  magnitude  of  the  product 
of  the  zeroth  and  +1  order  beams  (in  the  expression  for  the  H function) 
increases  as  w/v  increases  simply  because  the  deflected  beam's  axis 
moves  closer  to  the  detector.  In  fact,  the  product  is  maximum  when 
the  deflection  is  such  that 


However,  once  the  1 order  beam  is  deflected  out  of  the  region  for 
which  the  index  of  refraction  spectrum  is  constant,  the  scintillation 
spectrum  will  decrease  (since  $n  decreases  beyond  this  point).  There- 
fore the  scintillation  spectrum  will  increase  up  to  the  point  at 
which 


Lb.L  1.077 

r . - k i0 


U) 


1.077 

~LT“ 


and  will  decrease  beyond  this  point. 
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increases  as  u>/v  increases  simply  because  the  deflected  beam's  axis 
moves  closer  to  the  detector.  In  fact,  the  product  is  maximum  when 
the  deflection  is  such  that 


However,  once  the  1 order  beam  is  deflected  out  of  the  region  for 
which  the  index  of  refraction  spectrum  is  constant,  the  scintillation 
spectrum  will  decrease  (since  4>n  decreases  beyond  this  point).  There- 
fore the  scintillation  spectrum  will  increase  up  to  the  point  at 
which 


L u _ L 1.077 
¥ v " k Lq 


1 .077 
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When  the  detector  is  at  the  90  degree  location  (r0  = (0,yo)) 
the  < variable  of  integration  sweeps  the  deflected  beam  through  all 
possible  orientations  with  respect  to  the  detector.  Hence  the 
scintillation  spectrum  will  remain  constant  up  to  the  frequency 

a)  1.077 


and  decrease  thereafter. 

The  purpose  of  the  preceding  discussion  was  to  provide  a physical 
interpretation  of  the  effects  observed  when  the  detector  is  off-axis 
to  the  beam.  Of  course  if  the  detector  is  axially  located  we  should 
expect  to  see  none  of  the  behavior  which  has  just  been  described. 
Qualitatively,  the  results  of  this  discussion  and  of  the  earlier 
analytic  evaluation  are  identical,  thus  lending  credence  to  those 
results. 

Asymptotic  Behavior  of  S2  (o> ) 

The  S£  spectrum  of  Equation  (126c)  is  simply  the  spectrum  of  a 
gaussian  beam  with  an  axially  located  detector.  As  such,  it  de- 
termines the  behavior  of  the  spectrum  of  Equation  (126a)  for  an  axial 
detector  as  well  as  for  the  limiting  plane  and  spherical  wave  cases. 

We  shall  be  concerned  with  two  basic  frequency  ranges;  high  frequencies, 
such  that 


*2-  » 2tt, 
“O 


and  low  frequencies,  such  that 


2-  « 2tt. 
wo 

In  addition,  within  the  high  frequency  region,  we  are  interested  in 
two  sub-regions: 

N(u)/u>  )^ 

5 *-  > 1 (sub-region  A) 

2,[HN2(l+e)2] 


N(oj/u)  ) 

y 5—  < 1 (sub-region  B). 

2tt[1+N  (1+er] 
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Within  the  high  frequency  region,  i.e.  frequencies  such  that 


1 


— >>2*  , (153) 

wo 

the  sin2  function  is  oscillating  so  rapidly  with  respect  to  the  other 
terms  in  the  integrand  of  Eq.  (122c)  that  we  may  set 

sin2(- • • ) % 1 . (154) 


The  S2  spectrum  is  then  given  by 


S2(w)  = 4 die*  [ki2  + (to/w0)2] 

0 

N[k' 2+(u/w0)2]  2 

x f1  du  e 2ir[l+N2(l+e)Z]  U (155) 

J0 

where  we  have  made  use  of  the  change  of  variable  u=l-z.  For  very 
high  frequencies  (high  frequency  sub-region  A)  the  u integrand  is 
small  long  before  the  upper  limit  of  u=l  is  approached-  Therefore 
the  upper  limit  may  be  extended  to  infinity  and  the  u integration 
performed  in  closed  form  to  yield 


N/2ir  “H/2  « 2-i"7/3 

So(co)  = 2 /T  5 5-  die'Dc'*  + (u>/o>0)  J - (156) 

L1+N2(1+£)^J  jo 


The  change  of  variable 
x = k,2/(u/w0)2 

puts  Eq.  (156)  in  terms  of  a Beta  function  [73]  integral  which  gives 
for  the  S2  spectrum 
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r , ^ r(l  1/6)  r N/2tt 
s2(“»-  -rW  T^Tjl 


-1/2  , -11/3 


(157) 


For  the  case  in  which  the  upper  limit  on  the  path  integration 
cannot  be  extended  (high  frequency  sub-region  B) 


S2(u)  = 2 ^ ^du  e'bu  d<  < 1//2  (k+1) 


-11/6  e’bu  < 


where 


(N/2ir)(«/«o  )c 

b ? — S- 

1 + r(l+0 

and  we  have  used  the  change  of  variable 


k = k 1 2/  (to/ 


The  k integration  is  of  the  form  [74] 


dt  t2v-]  (t+l)2^"1  e'Pt 

J n 


= r(2v)  p"y"v  eP/2  W + .(p)  • 

' r y-v,y+v-Js 


(158) 


(160) 


where  W is  Whittaker's  function  [75], [76]  which  can  in  turn  be  expressed 


Wv-v,y-v-%(P)  = 4K-Sr^  PP+V  e P/2  1^1  (2v;2y+2v;p) 


+ rl2jM+2v^lJ_  pl-M-v  e-p/2  ^F]  (l-2y  ;2-2y-2v;p) , 


(161) 


S 


where  7 is  Kummer's  confluent  hypergeometric  function  [76]  defined 
as 


, x " r(x+n)  r(y)  zn 

(x;y;z)  = L rTxTrW  nr 


With  these  formulae  Eq.  (158)  is 


(162) 


s2<“>  ■ 2 ^ (^)'8/3  /’*■  e'bu2  M --W' 


2 r(-4/3)^"8/3  |3du  e-bu^(bu2)4/3  -j F, (li;Z.;bu2) . (163) 


Employing  Eq.  (162)  and  integrating  term  by  term  yields  finally  [ 771 [ 78] 

V8/3  l A r(4/3)  ; d/2)n  [?%)] 

S2(to)  ~ U0J  } rtlW  l <-1/3)n  nT(n+l/2) 


x 1F1 


+ r(-4/3)  l 
n=0 


' 1 3 «r/ii  \z\ 

n+?;n+2;  - 2^v*r J J 

fer/«!fln+4/3  2 

177317  mTn-Mi/e) 17  6 ’n^s’MvJ  ’ 


(164) 


where  (a)n  is  Pochhammer's  symbol  defined  as 


and  ar  is  defined  as 


^ l+N2(l+e)2 


(165) 


(166) 
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By  integrating  the  hypergeometric  functions  in  Equation  (163)  term  by 
term  we  have  made  the  tacit  assumption  that  these  series  are  uni- 
formly convergent  [79].  That  this  is  so  may  be  shown  by  calculating 
the  radius  of  convergence  [80]  of 


,F.(a;b;z)  = l C zn, 
11  n=0  n 

where 


(167a) 


(167b) 


The  radius  of  convergence,  R,  is  given  by 


lim 

OX” 


1 

R . 


For  the  confluent  hypergeometric  function 


i = lim 

n-xo 


(a+n) 

fb+n) (1+n) 


= 0, 


(168) 


(169) 


i .e. , 


R . (170) 

Therefore  the  (confluent)  hypergeometric  series  is  uniformly  con- 
vergent [81]  for  all  values  of  the  argument  z (which  is  in  general, 
complex),  and  we  are  justified  in  integrating  the  series  term  by 
term. 

If  we  retain  only  the  first  few  terms  of  the  series  in  Equation 
(164)  the  behavior  of  S2  within  the  high  frequency  region  may  be  more 
closely  examined; 


S2(w) 


2 JTTr  (4/3) 
(11/6) 


r 


r(-4/3)r(l 1/6) 
yir(4/3) 
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The  dominant  behavior  in  this  frequency  range  is  seen  to  be  the  -8/3 
power  of  frequency. 

The  low  frequency  behavior  of  S2  is  now  calculated  by  retract- 
ing the  approximation  of  Equation  (154),  letting 


sin2x  = \~\  Re(el2x}  , 


(172) 


and  employing  techniques  and  formulae  similar  to  the  preceding  [74,82, 
83].  This  approach  results  in  the  formula 

c / \ or- .-an/l  r (4/3)  “ (]/2)n 

S2(o))-2^TA  <2  F(ii/6)  lQ  (-l/3)n  n! 

(n+  2';n+^;-cxr,B2/2iT)j  p(_4/3)  “ ( 1 1 /6 )n  (a  A2/2iT)n+4/^3 


iB2\m  2fi (-(m+n) ,m+n+l ;m+n+2;-ia) 


F,  (-(m+n+  o-),m+n+  o-;m+n+ 


m+n+7/3 


I 


where 


“ ~ l-iN(l+e) 


N 


ar  = Re{a)  - 0 0 

- ■ Mt /i  . 


l+N6(l+e) 


A = 


tef  • W 


1/2 


B = oj/w. 


(173b) 

(173c) 

(173d) 

(173e) 


and  2^1  is  Gauss'  hypergeometric  function  [84]  defined  as 

» (w)  (x)  n 

2^1  (WtX;y;z)  = l ~JZ) fr 

n=0  vy'n  n* 


(174) 


While  we  have  made  no  approximations  in  deriving  Eq.  (173a),  its 
usefulness  is  limited  to  low  frequencies  because  the  series  do  not  con- 
verge uniformly  for  N on  the  order  of  unity  or  for  high  frequencies. 

If  only  the  first  few  terms  of  the  equation  are  retained,  we  obtain 


S2(u>)  = 


2^  r( 


(4/3)  |~27  n \4/3  r(l  1/6)  /F1P(N) 
r(l  1/6)  2Z\ft)  r (7/6)  Vf1s(N), 


+ F?(N)(A"8/3  B4  + 3A"2/3  B2  - |^A4/3) 


(175a) 


where  the  function 


F1P(N)  = aj/3  + jj-  1 " (“i  + & ar^ 
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(a2  - + 2/3  araj) 


(175b) 


gives  the  correct  plane  wave  D.C.  level  as  N-+O,  and  the  function 
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(175c) 


F,S(N) 


a4/3  tJ?r(4/3) 
r 3r(ll/b)2  ' 


goes  to  the  spherical  wave  D.C.  level  as  . Also  in  Equations  (175) 
we  have  defined 


and 


a..  = Im{a} 


(1 75d) 


(1 75e) 


To  within  a multiplicative  constant  (due  to  the  normalization  used 
herein)  Equation  (175a)  can  be  shownto  yield  the  limiting  plane  and 
spherical  wave  (N-*0,N-*>°)  behavior  as  quoted  by  Reinhardt  and  Collins 
[31]. 


As  an  indication  of  its  plane  and  spherical  wave  limits  we  have 
displayed  in  Figure  11  the  D.C.  level  of  S2  normalized  to  the  spherical 
wave  D.C.  level,  as  a function  of  the  Fresnel  number  N.  Although 
Figure  11  indicates  a lack  of  agreement  between  expressions  (175b)  and 
(175c)  within  the  region  NM , the  figure  is  indicative  of  the  con- 
ditions under  which  the  gaussian  beam  may  (for  the  purposes  of  ob- 
serving scintillation  spectra)  be  considered  to  be  in  the  near  or  far 
fields. 

Final  results  of  the  asymptotic  evaluation  of  the  spectrum  S2  are 
presented  in  Figure  12a.  The  breakpoint  between  the  zero  slope  D.C. 
level  and  the  -8/3  slope,  which  is  commonly  referred  to  as  the  Fresnel 
breakpoint,  will  move  slightly  as  the  Fresnel  number  N of  the  gaussian 
beam  varies  between  zero  and  infinity.  However,  owing  to  the  steepness 
of  the  cutoff  (-8/3),  the  Fresnel  breakpoint  will  always  (for  weak  tur- 
bulence) be  on  the  order  of  2ir.  The  Fresnel  breakpoint  can  be  inter- 
preted as  denoting  the  temporal  coherence  length  (in  the  presence  of 
the  turbulence)  for  the  source  employed.  At  this  breakpoint  (at  least 
for  plane  and  spherical  waves) 

U>  „ o 

— 'v  2it 

u>_ 


or 
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FRESNEL 


Figure  11 —Variation  with  Fresnel  number  of  the  S2 
spectrum  D.C.  level. 
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Figure  12a--Asyinptotes  for  S£  spectrum. 

The  temporal  coherence  length  is  then  simply 
tQ  = 1/f  - JIl/v. 

An  expression  for  the  breakpoint  between  the  -8/3  and  -11/3 
asymptotes  is  obtained  by  setting 


This  breakpoint  frequency  is  minimum  for  the  collimated  laser  beam 
NM,  and  goes  to  infinity  in  both  limits  as  N-K)  and  From  Equa- 

tions (1 1 5d)  and  (149a)  this  point  is  seen  to  be  at 


w IZ  3 
7 r N2*  a 


r(jl/6) 

r (4/3 ) 
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W(L) 

JIl 


or  at 


. V_  2W(L) 

JTl  JIl 


where  W(L)  is  the  receiver  plane  beam  spotsize.  The  frequency  of  the 
break  point  is  increased  (decreased)  as  the  number  of  Fresnel  zone 
size  patterns  contained  within  an  area  the  size  of  the  receiver  plane 
beam  spotsize  increases  (decreases).  Thus  the  term 

2W(L)//aT 

constitutes  a beam  size  factor. 

Figure  12b  presents  numerical  values  for  the  asymptotes  when 


N = 2 x 10“ 


(177) 


LQ  = 10  meters 


and 


JJL  = 2.516x10"^  meters, 


plotted  from  Equations  (157),  (171),  and  (175). 
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As  a final  result  of  the  preceding  asymptotic  evaluation  we 
present  in  Figure  13  a composite  of  the  S]  and  S2  spectra.  The 
slopes  and  breakpoints  have  been  given  earlier.  In  Figure  13  the 
departure  from  a flat  spectrum  for  frequencies  below  the  Fresnel 
breakpoint  is  a result  of  the  detector  being  off-axis.  Previous 
authors  [37], [39]  who  always  assumed  an  axially  located  detector  pre- 
dicted no  such  behavior.  Other  authors  [85], [86]  assumed  that  low 
spatial  frequency  turbulence  or  system  noise  (for  the  case  in 
which  the  beam  is  being  pointed  by  a servo-system)  gives  rise  to  beam 
wander,  which  is  equivalent  to  saying  that  r0  is  a random  variable. 
However  these  workers  were  not  concerned  with  the  resulting  effect 
upon  the  temporal  scintillation  spectra. 

This  concludes  the  discussion  of  the  asymptotic  evaluation  of 
the  gaussian  beam  scintillation  spectrum. 


Figure  13— Asymptotes  of  gaussian  beam  spectrum. 


D.  Numerical  Evaluation 

This  section  is  concerned  with  a numerical  evaluation  of  Equation 
(126)  for  the  temporal  scintillation  spectrum.  The  computer  programs 
were  written  in  Fortran  IV  and  employed  Gaussian  quadrature  inte- 
gration techniques.  Listings  of  these  programs  plus  brief  descrip- 
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tions  are  contained  in  Appendix  G.  In  addition  to  the  numerically 
calculated  spectra,  we  also  present  computer  generated  plots  of  the 
normalized  differential  path  contribution  to  the  spectra,  defined 


as 


S» 


-1 


r(o,o,L)o;-i 


d 

Si 


a,oSD(a)) 

~TT — 


Lr(0,0,L)o‘J 


(178) 


Presented  in  Figure  14  are  the  computer  calculated  spectra  of  a 
clean  gaussian  beam  for  various  cases.  Note  the  similarity  between 
these  curves  and  the  curves  in  Figures  12a  and  13.  The  decrease  in 
the  magnitude  of  the  spectra  for  the  off-axis  cases  is  due  simply 
to  the  decreased  intensity  of  the  beam  at  these  points.  This  figure 
shows  the  spectrum  resulting  from  a laser  beam  with  Fresnel  number 
N 'v  2 x 103.  Curves  resulting  from  axial  as  well  as  off-axis  point 
detectors  are  included.  Note  that  for  r<)  = (x ' ,0 ) the  low  frequency 
peak  is  at  0 


u/«0  (1 .077/Lq) 


Now  consider  the  differential  path  contribution  to  the  spectrum. 

To  put  things  in  perspective,  Figure  15  shows  the  computer  generated 
differential  path  contributions  for  plane  and  spherical  waves  for 
w/uo  < 2tt . For  the  plane  wave,  the  differential  path  contribution  is 
peaked  at  the  transmitter.  In  contrast,  the  spherical  wave  dif- 
ferential path  contribution  is  peaked  at  mid-path.  Figure  16  is  the 
same  except  for  u>/u>0  > 2tt.  In  this  frequency  range,  both  the  plane 
and  spherical  wave  differential  path  contributions  are  oscillatory. 
These  figures  are  to  be  compared  with  Figure  17  which  shows  the 
gaussian  beam  differential  path  contributions  for  o)/u>0  < 2tt  and 
w/wo  > 2it.  For  this  case  the  Fresnel  number  was  on  the  order  of  2 x 
103.  From  the  results  of  the  analysis  within  this  chapter  it  is  seen 
that  the  gaussian  beam  should  behave  as  a spherical  wave,  yet  when  the 
detector  is  off  axis,  the  low  frequency  differential  path  contribution 
appears  to  be  that  of  a plane  wave.  That  is  the  path  contribution  is 
more  heavily  weighted  towards  the  transmitter.  This  phenomenon  has 
important  implications  in  situations  where  laser  beams  are  propagated 
between  aircraft.  For  an  off-axis  detector  the  scintillation  spectrum 
will  become  sensitive  to  the  outer  scale  of  any  turbulence  localized 
in  the  vicinity  of  the  transmitting  aircraft's  fuselage.  For  suffi- 
ciently short  outer  scales  (typically  £ 1 meter)  [87]  the  Fresnel  break 
point  at  ui/o)  2tt  may  easily  be  obscured. 
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Figure  14— Typical  gaussian  beam  spectra  for  axial 
and  off-axis  detectors. 
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Figure  15— Low  frequency  differential  path  contribution 
for  plane  and  spherical  waves. 
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Figure  16— High  frequency  differential  path  contribution 
for  plane  and  spherical  waves. 
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Figure  17— Differential  path  contribution  for  gauss i an 
beam  with  off-axis  detector. 


This  concludes  the  analysis  for  a clean  gaussian  beam  with  a 
point  receiver.  In  the  next  section  we  will  extend  the  results  of 
this  and  the  previous  section  to  include  the  effects  of  finite  de- 
tector apertures  and  a finite  inner  scale. 

E.  Additional  Topics 

We  now  wish  to  extend  the  results  of  Section  C to  account  for 
two  phenomena  not  previously  mentioned,  namely  the  effects  upon  the 
scintillation  spectrum  of  finite  detector  apertures  and  a finite 
inner  scale.  The  former  effect  is  of  interest  because  all  physically 
realizable  detectors  have  a finite  aperture.  The  assumption  in 
Section  D of  a point  detector  was  merely  a device  used  to  simplify 
the  mathematics.  However,  we  will  demonstrate  that,  under  certain 
conditions,  this  is  a reasonable  approximation  to  reality.  Inclusion 
of  a finite  inner  scale  in  the  model  of  the  index  field  reflects  the 
fact  that  the  kinetic  energy  of  very  small  scale  index  inhomogeneities 
is  dissipated  in  the  form  of  heat  rather  than  being  passed  on  to 
ever-decreasing  scale  sizes.  From  the  discussion  of  phase  gratings 
in  Chapter  II,  it  is  seen  that  the  existence  of  a finite  inner  scale 
requires  that  (for  weak  turbulence)  the  angle  of  the  bending  ex- 
perienced by  light  rays  passing  through  the  turbulent  atmosphere  be 
bounded  by 

“max  = X/V 

where  A is  the  wavelength  and  Jt0  is  the  inner  scale  of  turbulence. 


In  Eq.  (123)  of  Section  C we  could  have  easily  included  an  inner 
scale  cut-off  in  the  expression  for  the  power  spectrum  of  the  index 
fluctuations  [71]; 

-(k/k  )2 

*n(|e1»lc2,,c3)  = 0,033  Cn^2  + O-077/Lo)2]"11/6  e m , (179) 


where 


2 2 2 , 2 

k = <i  + <2  + <3 


Km*0  - [0. 033ttT (5/3 ) ]~3/4  = 5.92 


and  Ji0  is  the  inner  scale  of  turbulence. 


91 


Inclusion  of  this  inner  scale  dependence  results  in  an  additional 
term  in  the  integrand  of  Eq.  (125)  of  the  form 


exp 


-[k,2+Uu»0)2] 

(/AUJ2 


(180) 


Inspection  of  this  term  shows  that  the  inner  scale  will  produce  a 
gaussian  roll-off  of  the  scintillation  spectrum  in  the  vicinity  of 


y-%c92^k  . (181) 

“o  o 


-2  -3 

For  a Fresnel  length  on  the  order  of  the  inner  scale  Uo^lO  -10  m) 
this  effect  could  easily  dominate  the  high  frequency  behavior  of  the 
scintillation  spectrum.  If  for  example  the  range,  L,  is  one  kilo- 
meter, the  source  is  a helium-neon  laser  (A=.6328xl0*6m)  and  the 
inner  scale,  s,0 » Is  one  millimeter  then  the  scintillation  spectrum 
will  exhibit  an  inner  scale  effect  for 

— > 150  » 2tt. 
uo  * 

Or  for  a velocity,  v,  of  150  meters/sec.,  inner  scale  effects  would 
begin  at 

f * v/A0  = 150  KHz. 


Finite  aperture  effects  are  determined  by  defining  the  receiver 
aperture  as 


W(r)  = 


r 

i 

< 

o 

s. 


I r | < D/2 
|r|  > D/2 


(182) 


where  D is  the  diameter  of  the  photosensitive  area  of  the  detector. 
For  a spherical  wave  (chosen  because  the  integrations  are  easily 
performed  in  closed  form)  the  H function  is  given  by  (see  Equation 
(89b) 
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h(z.k'.W)  = |E0|2  e 2k 


|dr  W(r)  eiz< 


" IE0!  e 


D/2os  „ f2"  d»  e’zK  Dcos^-e^  (183) 

0 J0 


or  finally 


H(z,k\W)  = 


\2  ? 2J1  (z  2 K0 

)lE°l 


- i z(l-z)L<“ 
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(184) 


The  scintillation  spectrum,  from  Eq.  (92),  is  then 


SDU)  a l\z  d<  [,c2+(aJ/to_)2+(l  .077/L’ )23 

H Jo  Jo 


-11/6 


„.2  fz(l-z)  r 2./  , C2Jl(z  2 J<2+(“/“  0)2) 

sin  1-V-  l*  +(u/%)  D * I 2~.  f 72 

L i 7 z V +WwoJ 


(185) 


Within  the  high  frequency  region  we  can  approximate  the  sin2  function 
by  1/2,  extend  the  z integration  to  infinity[88] , and  perform  the  k 
integration  to  obtain 


V-1 » j fc) 


-11/3 


(186) 


Similar  results  are  to  be  expected  for  plane  or  gaussian  beam  waves. 
The  onset  of  this  so  called  aperture  averaging  phenomenon  may  be  de- 
termined by  equating  the  formulae  for  the  -8/3  asymptote  (derived 
under  the  assumption  of  a point  detector)  and  the  -11/3  asymptote 
for  a finite  aperture. 
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This  aperture  averaging  region  can  also  be  discerned  by  physical 
arguments.  At  the  receiver  we  observe  a series  of  diffraction  pat- 
terns of  various  sizes,  I,  flowing  in  the  direction  of  the  wind  at 
velocity  v.  The  associated  temporal  frequency  at  the  output  of  the 
detector  will  therefore  be 

f = v/t 


or 


O)  _ 2tt 

w0  l/J\ L 


(187) 


However  for  an  aperture,  0,  larger  than  the  dimension  l , the  detector 
wil^  integrate  over  these  diffraction  patterns  [62]  thus  attenuating 
the  effects  of  this  scale  size.  Aperture  averaging  would  then  occur 
for  frequencies  such  that 


U) 


> 


2ir 

d/JTl 


(188) 


Therefore  aperture  averaging  should  be  of  little  concern  if  we  require 
that  the  aperture  diameter  be  much  smaller  than  the  Fresnel  zone; 

D <<  J XL  . 


As  an  example,  for  a one  kilometer  range,  a HeNe  source,  and  a one 
centimeter  receiver  aperture,  one  would  expect  to  observe  aperture 
averaging  for 


Alternatively,  for  a velocity  of  150  meters/sec.,  aperture  averaging 
would  occur  for  frequencies 

f > v/D  = 15  KHz. 


On  the  basis  of  arguments  presented  in  this  section  it  is  seen 
that  if  we  require 


&0»D  « JTl  , 


(189) 
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then  the  inner  scale  effects  and  aperture  averaging  effects  should  be 
of  little  concern  in  the  determination  of  the  scintillation  spectrum. 
This  concludes  the  analysis  of  the  temporal  scintillation  spectrum 
of  a clean  gaussian  beam. 

F.  Summary  and  Conclusions 

Within  this  chapter  we  have  presented  an  evaluation  of  the 
temporal  scintillation  spectrum  of  a clean  gaussian  beam  fur  the 
case  in  which  the  detector  is  off-axis  to  the  beam.  The  off-axis 
detector  problem  has  not  been  previously  treated  in  the  literature. 

The  general  formula  for  the  spectrum  derived  in  Chapter  II,  was 
employed  in  conjunction  with  a particular  model  of  the  index  of 
refraction  spatial  spectrum,  plus  the  assumption  of  an  off-axis 
point  detector,  in  deriving  an  expression  which  was  then  evaluated 
asymptotically  and  numerically.  As  a result  of  the  analysis  it  was 
shown  that  the  fact  that  the  detector  was  off  axis  produces  some  in- 
teresting outer  scale-dependent  behavior  for  frequencies  below  the 
Fresnel  breakpoint.  In  addition,  this  low  frequency  behavior  was 
shown  to  be  dependent  upon  the  location  of  the  detector  with  respect 
to  the  beam  axis  and  wind  direction. 

„ ^Section  B was  concerned  with  setting-up  the  problem.  The  von- 
Karman  model  of  the  index  spectrum  together  with  a gaussian  beam  source 
field  were  used  in  the  expression  for  the  spectrum  derived  in  Chapter 
II.  The  assumption  of  a point  detector  enabled  performance  of  the 
integration  for  the  H function  of  Chapter  II.  The  final  expression 
for  the  scintillation  spectrum  was  demonstrated  to  reduce  to  that 
which  was  derived  by  other  authors  when  the  detector  was  axially 
located,  as  well  as  to  the  expressions  for  the  plane  and  spherical 
wave  scintillation  spectra  when  the  beam  Fresnel  number  was  allowed 
to  go  to  zero  and  infinity  respectively. 

In  Section  C the  expression  for  the  spectrum  was  evaluated.  The 
first  portion  of  the  section  was  devoted  to  developing  analytic  ex- 
pressions for  the  spectrum  of  a gaussian  beam  with  arbitrary  Fresnel 
number.  Expressions  for  the  spectrum  in  terms  of  hypergeometric 
functions  were  developed  and  asymptotic  behavior  (for  high  and  low 
frequencies)  was  determined  by  retaining  only  the  first  few  terms 
of  the  series  expressions  for  these  functions.  Under  certain  con- 
ditions it  was  found  possible  to  determine  closed  form  expressions 
for  the  spectrum  asymptotes. 

Section  D was  devoted  to  a computer  evaluation  of  the  spectrum 
for  a particular  beam  Fresnel  number,  and  a variety  of  detector 
locations.  In  addition,  differential  path  contribution  plots  were 
presented  for  plane  and  spherical  waves,  and  for  gaussian  beams 
with  off-axis  point  detectors. 
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As  a result  of  this  section  it  was  shown  that  for  the  detector 
and  beam  axis  along  a line  collinear  with  the  wind  direction,  the 
scintillation  spectrum  displayed  a peak  at  io/io0  = /J(l  .077/Lq)  with 
a +2  power  law  for  frequencies  below  and  a -2/3  power  law  for  fre- 
quencies above  this  point.  For  the  detector  and  beam  axis  along  a 
line  perpendicular  to  the  wind  direction  the  spectrum  showed  a break- 
point at  u)/u>  =1 .077/L'  with  a -2/3  power  law  for  frequencies  above 
this  point. 

The  differential  path  contribution  plots  showed  that  despite 
the  fact  that  a gaussian  beam  may  be  well  within  the  Fraunhofer  zone 
(N  >1),  the  off-axis  detector  causes  the  path  contribution  to  be 
peaked  at  the  transmitter  (as  it  is  for  a plane  wave)  rather  than 
at  midpath  (as  it  is  for  a spherical  wave).  This  effect  has  important 
implications  (as  will  be  demonstrated  in  the  next  chapter)  for  cases 
in  which  the  turbulence  strength  is  non-uniform  along  the  propagation 
path. 

In  Section  E the  analysis  was  extended  to  describe  finite  inner 
scale  and  finite  aperture  effects.  Through  mathematical  and  physical 
arguments  it  was  demonstrated  that  if  the  inner  scale  and  the  aperture 
diameter  are  much  smaller  than  the  plane  wave  Fresnel  zone  (/IT)  then 
their  effects  can  be  ignored. 

Throughout  this  chapter  the  off-axis  distance  and  orientation  of 
the  detector  have  been  assumed  to  be  deterministic.  For  the  case  in 
which  the  laser  beam  is  being  steered  by  a servo  system  however,  these 
variables  could  easily  be  considered  random.  By  assuming  appropriate 
statistical  distributions  for  these  variables  it  would  then  be  possible 
to  determine  the  effects  upon  the  scintillation  spectrum  of  servo 
tracking  error  or  noise. 


96 


CHAPTER  IV 

THE  DIRTY  GAUSSIAN  BEAM 


A.  Introduction 

The  purpose  of  this  chapter  is  to  develop  a model  of  a dirty, 
i.e.,  spatially  corrupted  laser  beam,  and  to  determine  the  effect 
upon  its  scintillation  spectrum  of  the  laser  beam  being  dirty.  This 
objective  is  pertinent  because  laser  beams  sometimes  do  not  have  the 
perfect  spatial  properties  which  most  analytical  studies  presume. 

Even  if  the  beam  does  have  good  spatial  quality,  many  situations 
require  it  to  be  passed  through  a window,  for  example  on  an  aircraft. 
This  unavoidably  introduces  a certain  amount  of  deterministic  spatial 
noise  into  the  phase  and/or  amplitude  of  the  beam. 

The  configuration  chosen  for  analysis  within  this  chapter  is  that 
of  a clean  gaussian  beam  shining  through  a window  whose  transmission 
function  is  assumed  to  be  unity  except  within  a small  axially  located 
spot  where  the  phase  of  the  incident  beam  is  shifted  slightly.  This 
model  is  appropriate  for  the  situation  in  which  the  laser  beam  is 
transmitted  through  a window  on  which  there  is  an  oil  spot  or  an 
imperfection. 

Instead  of  calculating  the  actual  receiver  plane  fields  due  to 
such  a configuration,  (as  required  by  the  general  formula  for  the 
scintillation  spectrum  which  was  derived  in  Chapter  II)  the  field 
immediately  past  the  window  is  expanded  in  a two  dimensional  series  of 
functions  involving  Gauss-Hermite  polynomials.  The  particular  formula- 
tion used  enables  the  receiver  plane  fields  to  be  written  down  immedi- 
ately without  having  to  resort  to  performing  the  diffraction  integrals. 
Since  the  size  and  phase  of  the  spot  are  chosen  a priori,  and  we  wish 
to  use  as  few  terms  of  the  polynomial  expansion  as  possible,  an  optimi- 
zation procedure  for  approximating  the  fields  behind  the  window  is 
presented. 

The  integral  expression  for  the  scintillation  spectrum  resulting 
from  our  formulation  of  the  dirty  beam  problem  is  evaluated  numerically 
and  the  results  presented  as  a series  of  plots.  As  a result  of  the 
analysis  the  dirty  beam's  scintillation  spectrum  is  shown  to  be 
sensitive  to  turbulence  localized  in  the  vicinity  of  the  transmitter. 
This  sensitivity  gives  rise  to  low  frequency  behavior  similar  to  that 
demonstrated  in  Chapter  III. 

In  Section  B the  dirty  beam  representation  is  formulated.  The 
dirty  beam  is  modeled  as  a clean  gaussian  beam  which  has  been  trans- 
mitted through  a window  having  a unity  transmission  function  except 


97 


for  a small  axially  located  spot  where  the  transmission  function  is 
complex.  The  fields  to  the  immediate  right  of  the  window  are  expanded 
in  a complete  orthonormal  free-space  eigenmode  set  [90,91].  These 
functions  which  closely  describe  the  modes  of  an  optical  resonator 
(i.e.,  a spherical  mirror  laser)  are  exact  eigenfunctions  (within  the 
paraxial  approximation)  of  the  free-space  wave  operator.  Therefore  an 
arbitrary  field  may  be  expanded  in  some  plane  in  terms  of  these  func- 
tions. The  functions  then  can  propagate  independently  to  any  other 
plane  and  be  reassembled  with  the  known  weights  to  determine  the 
actual  field.  This  property  of  the  Gauss-Hermite  polynomial  expansion 
is  employed  in  the  evaluation  of  the  general  expression  developed  in 
Chapter  II. 

By  use  of  Babinet's  principle  [92]  the  dirty  spot  is  essentially 
expanded  in  terms  of  a truncated  series  of  laser  modes.  Since  this 
series  is  truncated,  we  present  a procedure  for  optimizing  the  approxi- 
mation. 

A particular  case  of  the  dirty  beam  problem  is  treated  in  Section  C. 
The  required  integrations  are  carried  out  numerically  for  a variety  of 
detector  locations  and  turbulence  strength  distributions  and  the 
results  presented  in  the  form  of  several  plots  of  spectra. 

Section  D contains  a summary  and  conclusions. 


B.  Development  of  Model 


This  section  is  concerned  with  the  formulation  of  a model  of  the 
dirty  beam  and  an  optimization  procedure  for  approximating  it  in  terms 
of  eigenmodes  of  the  free  space  wave  operator. 

We  shall  first  demonstrate  the  use  of  this  particular  eigenmode 
expansion  as  a means  of  expressing  an  arbitrary  field.  The  specifi- 
cation of  a general  field  in  terms  of  the  Gaussian-Hermite  polynomial 
functions  is  written  as  [91] 


£(x,y,z)  = I l fm(x,z+z0)  fn(y,z+z0) 
m=0  n=0 


(190-a) 


where  the  f functions  are  defined  as 


fm(x.z+2o)  - (Hr)  * n-iN(z*z„)/Lf/2 


[l+iN(z+zn)/L]"(m+1)/2  H (, 


-a(z+z0)x‘ 


, (190-fc) 
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a(z+ZQ)  = 


kN/2L 


ar(z+z0)  = Re{a(z+z0) } 


(190-c) 


Hm  is  the  order  Hermite  polynomial  [93],  N is  the  beam  Fresnel 
number  defined  in  Equation  (115b),  and  the  beam  waist  is  located  at 
z=-z0.  Making  use  of  the  orthonormal ity  [94]  of  the  mode  functions 


[ dx  fK(x,z+z0)  f£(x,z+z0)  = 6KL 

J -oo 

gives  the  expansion  coefficients 

oo 

' 

= dxdy  E(x,y,z)  f*(x,z+z0)  ffj(y,z+z0) 
J 


(191) 


(192) 


The  field  of  a spherical  mirror  laser  operating  in  the  lowest  order 
transverse  spatial  mode  nominally  consists  of  only  the  zero-zero 
eigenmode.  If  we  arbitrarily  assume 


JJI 


; m,n=0 


0 ; m,n^0 


then  the  field  is  simply 


(193) 


E(x,y,z)  = ^N/2L-  e 


kN/2L  , 2.  2x 
l+iN(z+z0)"/L  (x  +y  } 


(194) 


With  the  exception  of  the  trivial  phase  factor  exp[ik(z+z0)]. 
Equation  (194)  is  recognized  as  the  expression  for  the  field  used  in 
Chapter  III  (Equation  (115a)).  Throughout  the  remainder  of  this  chapter 
we  assume  that  the  unperturbed  laser  beam  is  described  by  Equation 
(194). 

We  now  proceed  in  the  development  of  a particular  model  for  the 
dirty  laser  beam. 
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With  reference  to  Figure  18,  we  define  a circular  region  at  the 
2-0  plane  (the  plane  of  the  dirty  window)  within  which  the  field  of 
Equation  (194)  is  corrupted. 


Figure  18— Illustration  of  phase  obstacle. 


To  the  immediate  right  of  the  z=0  plane  the  field  is 


E(x,y,0+)D6(x-xs ,y-ys) 


(195-a) 


where  the  window  transmission  function,  D6,  is  given  by 


D6(x-xs,y-y$) 


A eiw  ;J  (x-xs)2+(y-ys)2  £ 6/2 


1 ; elsewhere 


(195-b) 


6 is  the  diameter  of  the  circular  region  and  the  field  is  given  by 
Equation  (194).  Equivalently  by  use  of  Babinet's  principle  we  may 
express  the  field  within  this  plane  as 


E(x,y,0+)  - [l-D6(x-xs,y-ys)]E(x.y,0+) 


(196) 
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or  by  defining  a new  variable,  D^, 


E(x,y,0+)  - D^(x-xs,y-ys)E(x,y,0+) 


0 97) 


In  the  material  to  follow,  we  shall  define  the  complementary  field  as 


Ec(x,y,0+)  = D^(x-xs,y-ys)E(x,y,0+) 


098-a) 


where 


1-A  eiy  ; J (x-xs)2+(y-ys)2  <6/2 


; elsewhere 


(198-b) 


As  an  example  consider  the  case  in  which  A=0.  Then  the  comple- 
mentary field  is 


r . E(x,y,0+)  ;J  (x-x  )2+(y-ys)2  < 6/2 

EC(x,y»0+)  =]  5 

0 ; elsewhere 
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Within  the  region  to  the  immediate  right  of  the  disc,  the  complementary 
field  exactly  cancels  the  incident  field. 

The  eigenmode  expansion  will  now  be  applied  to  the  complementary 
field: 


Ec(x,y,z)  = l l Ymn  gm(x,z+z-i)gn(y.z+zi)  » (200-a) 

m=0  n=0 


where  the  mode  functions  are  given  by 


x [l+iM(z+z-| )/L] 


-(m+1 )/2 


-efz+z^x* 


(200-b) 


I 


B(z+z,)  - W2  L- 

I 1 + lM(z+2-|  )/L 


6r(z+z-|)  s Re{e(z+Z1)} 


M is  the  Fresnel  number  of  the  complementary  beam  defined  as 


(200-c) 


M = 


XL 


(200-d) 


and  W-j  is  the  spotsize  of  the  beam  waist  which  is  located  at  z=-Z] . 


In  Appendix  F we  approximate  the  complementary  field  by  the  zero- 
zero  term  of  its  eigenmode  expansion.  This  approximation  is  optimized 
by  minimizing  the  integral  square  difference  between  the  actual  comple- 
mentary field  and  its  truncated  series  approximation. 

Also  in  Appendix  F we  justify  the  claim  that  without  loss  of 
generality  we  can  assume  the  phase  obstacle  to  be  axially  located  with- 
in the  incident  (main)  beam.  This  is  subject  to  the  conditions  that 
the  detector  be  near  the  center  of  the  main  beam  (r0  « W(L))  and  that 
the  receiver  be  in  the  far  field  of  the  main  beam. 

As  shown  in  this  appendix  the  complementary  beam  mode  expansion 
coefficient,  y00*  9iven  by 


yoo 


W0W1 ( 1 “ A eip) 


o(z0)b*(z1) 
a(z0)+B*(z1 )' 


X HI 


-[a(z0)+e*(z1)](«/2) 

e 


(201 -a) 


and  the  complementary  beam  Fresnel  number  and  beam  waist  location  are 
given  by  the  solution  of  the  simultaneous  equations 

Im{a(z0)}  = Im{ 3 ( Z] ) } (201 -b) 


and 
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(201 -c) 


where 


x 


Sr(2l> 

w 


and 


y = ar(z0) 


(201 -d) 


(201-e) 


With  the  complementary  beam  thus  completely  specified,  we  may  proceed, 
in  the  manner  of  the  development  of  Chapter  III,  to  calculate  the  H 
function.  The  exception  is  that  the  expression  for  the  field  in 
Equation  (89-b)  is  replaced  by 

E(x,y,L)  - E^(x,y,L)  . (202-a) 

In  our  case  1=0,  and  the  total  field  assumes  the  form  of  the  difference 
of  two  gauss i an  beams: 


E(x,y,L)  - E^(x,y,L)  = 

-a(L+z0)(xV) 

V(L+zo)e 


-ed+z^Kx^y2) 

W^d+z-^e 


x 2WqW-| 


a(z0)6*(z1) 

a(z0)+3*(z1) 


(1-Aeip 


) f e-[a(z0)+B*(2l)](5/2)2| 


(202-b) 


This  completes  the  discussion  of  the  modeling  of  the  dirty  beam 
and  the  optimization  of  the  field  expansion.  These  results  will  be 
applied  in  the  next  section  to  the  calculation  of  the  scintillation 
spectrum  of  a particular  dirty  beam. 


C.  Evaluati on 

The  modeling  of  the  previous  section  was  applied  in  the  calcu- 
lation of  the  scintillation  spectrum  for  a particular  configuration; 
that  of  a laser  beam  with  a central  phase  perturbation.  The  spectrum 
resulting  from  this  field  configuration  was  calculated  for  three  tur- 
bulence models;  homogeneous  turbulence  with  a large  outer  scale, 
localized  turbulence  with  a small  outer  scale,  and  a combination  of 
these  two.  The  detector  was  assumed  to  be  off-axis. 
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We  evaluated  the  spectrum  resulting  from  the  situation  In  which  a 
laser  beam  with  waist  spotslze  (W0)  of  20.6  x 10"6  meter  was  directed 
through  a dirty  window  toward  a receiver  approximately  one  kilometer 
distant.  The  beam  waist  was  located  one  meter  Inside  the  window,  and 
the  beam  was  assumed  to  have  propagated  In  free  space  until  It  passed 
through  the  window  Into  the  turbulent  atmosphere.  The  phase  and  dia- 
meter of  the  spot  on  the  window  were  chosen  so  that  Its  receiver  plane 
diffraction  pattern  consisted  of  a dark  spot  approximately  one  meter  In 
diameter.  These  parameters  were  respectively  y»5.336  radians  and  «■ . 91 7 
x 10"3  meter.  The  values  of  these  variables  were  deduced  by  numerically 
calculating  the  receiver  plane  diffraction  patterns  of  obstacles  of 
various  phases  and  diameters.  Although  this  procedure  supplied  the 
receiver  plane  fields  as  required  by  the  general  formula  for  the  scin- 
tillation spectrum,  the  computations  were  lengthy  and  time  consuming. 
Therefore  we  resorted  to  the  approximation  of  the  complementary  field 
as  outlined  In  the  previous  section. 

The  routine  of  the  previous  section  gave,  for  the  complementary 
beam,  a waist  spotslze,  Wj,  of  3.146  x 10-4  meter  and  a beam  waist 
location,  -zi,  of  -.408  meter.  For  the  Helium-Neon  laser  (x».6328  x 
10’6  meter)  the  complementary  and  main  beam  Fresnel  numbers  were 
respectively  2.04  x 103  and  4.75  x 105,  and  the  receiver  plane  spot- 
sizes  were  W](L)  » .64  meter  and  W0(L)  ■ 9.79  meters.  The  parameters 
of  interest  are  Illustrated  In  Figure  19. 


Figure  19--Sketch  of  transmitter  plane  fields. 
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With  the  beam  waist  locations  and  spotsizes  for  the  main  and  com- 
plementary beams  determined,  the  fields  required  (see  Equation  (202b)) 
in  the  calculation  of  the  H function,  and  hence  the  spectrum,  were 
specified.  In  calculating  the  H_function  (Equation  (89b))  a point 
detector  at  a general  location,  r=r0,  was  assumed. 

A digital  computer  was  employed  in  performing  the  path  and  fre- 
quency integrations  of  Equation  (92).  For  the  path  integral  a single 
eight  or  32  point  Gaussian  quadrature  [95,96,97]  routine  was  used. 
Composite  integration  [98]  using  eight  point  Gaussian  quadrature 
routines  was  used  for  the  frequency  integral.  A listing  and  brief 
discussion  of  this  program  is  contained  in  Appendix  G. 

Figure  20  shows  the  spectrum  resulting  from  the  dirty  beam  for 
an  outer  scale,  L0,  of  100  meters  and  the  detector  positioned  at  r0  = 
(0,.88m).  This  spectrum  is  indistinguishable  from  the  spectrum  of  a 
clean  gaussian  beam  (see  Figure  14  in  Chapter  III).  Little  effect  of 
the  beam  being  dirty  was  to  be  expected  because  of  the  extreme  distance 
of  the  detector  from  the  diffraction  pattern  of  the  phase  obstacle 
(r0/W](L)=1.4),  and  because  the  detector  was  still  close  to  the  axis 
of  the  main  beam  (r0/W0(L)«l ). 

In  Figure  21  we  present  the  spectrum  resulting  from  the  same 
conditions  as  for  Figure  20  except  with  the  detector  located  at  r0= 
(0..44)  (r0/W-|(L)=.7).  The  low  frequency  behavior,  which  is  character- 
ized by  the  break  point  in  the  vicinity  of  u/u)0  = 10"4,  is  the  same  as 
that  shown  in  Chapter  III  for  the  off-axis  detector.  In  this  case 
however,  the  low  frequency  behavior  was  due  to  the  detector  being  off- 
axis  to  the  complementary  beam,  not  the  main  beam. 

The  spectrum  resulting  from  a different  turbulence  distribution 
is  presented  in  Figure  22.  In  this  case  the  outer  scale  was  one  meter 
and  the  detector  location  was  the  same  as  in  Figure  21,  however  the 
turbulence  structure  parameter  was  specified  as 


'°iZ1 


(203) 


where  z^  was  10  meters.  This  turbulence  strength  distribution  would 
be  appropriate  for  a situation  in  which  all  the  turbulence  was  confined 
to  a 10  meter  thick  region,  i.e.,  boundary  layer  at  the  transmitter. 

From  the  discussion  in  Chapter  III  of  the  spectrum  differential 
path  contribution  we  deduce  that  the  complementary  beam  (because  of 
the  off-axis  detector)  interacts  more  strongly  with  the  turbulence 
localized  near  the  transmitter  than  does  the  main  beam.  This  thought 
is  reinforced  by  the  order  of  magnitude  difference  between  the  low  fre- 
quency portion  of  the  spectrum  in  Figure  22  (due  primarily  to  the 
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Figure  20--Spectrum  for  detector  off-axis  by  a large  distance. 
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f 


NORMALIZED  FREQUENCY  uyw# 

Figure  21— Spectrum  for  detector  near  diffraction  pattern 


complementary  beam)  and  the  portion  of  the  spectrum  to  the  Immediate 
left  of  the  Fresnel  breakpoint  (due  to  the  main  beam). 


A more  realistic  model  of  the  turbulence  field  would  contain 
homogeneous  turbulence  as  well  as  a boundary  layer  at  the  transmitter. 
This  would  correspond  to  the  situation  In  which  the  laser  beam  was 
being  directed  between  two  aircraft.  In  fact,  to  be  totally  realistic, 
we  should  Include  a boundary  layer  at  the  receiver.  However,  the  dif- 
ferential path  contribution  plots  of  Chapter  III  have  shown  the  Insen- 
sitivity of  the  scintillation  spectrum  to  such  a turbulence  distribu- 
tion. Therefore  the  receiver  plane  boundary  layer  can  be  Ignored. 


Figure  23  shows  the  spectrum  resulting  from  a combination  of 
localized  and  distributed  turbulence  strengths.  For  this  case  the 
turbulence  field  model  used  was 


+ (1 . 077/Lq-j  )2]-11/6 


+ 0[k2  + (1 .077/L(!)2)2:-11/6 

4 


0 < z < 1 


(204) 


where  the  outer  scales,  Lq]  and  Lq?,  were  respectively  one  and  100 
meters,  the  e"1  point  of  the  localized  turbulence  was  ten  meters  and 
the  weights  a and  0 were  chosen  subject  to  the  restraint 


(205) 


The  restriction  of  Equation  (205)  was  Imposed  merely  to  facilitate 
direct  comparison  of  Figure  23  with  the  previous  results  of  this 
chapter  and  Chapter  III.  It  requires  the  path  integral  of  the  tur- 
bulence strength  to  be  a constant.  The  turbulence  model  of  Figure  22 
also  satisfies  Equation  (205).  For  the  spectrum  in  Figure  23  we  chose 
0s. 05.  Equation  (205)  then  gave  a*95,  for  a ratio  of  roughly  three 
orders  of  magnitude  between  the  turbulence  strength  in  the  boundary 
layer  and  that  In  the  distributed  turbulence.  Also  indicated  in 
Figure  23  are  the  two  outer  scale  and  two  Fresnel  breakpoints  which  are 
due  to  the  two  turbulence  strength  distributions.  This  figure  is 
Illustrative  of  the  Increased  sensitivity,  of  the  dirty  beam  scintil- 
lation spectrum,  to  an  Inhomogeneous  turbulence  distribution. 

Shown  In  Figure  24  Is  the  spectrum  resulting  from  a different 
combination  of  localized  and  distributed  turbulence.  For  this  case 
we  used 
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Figure  23--Spectr*um  for  combination  of  localized  and 
distributed  turbulence. 
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Figure  24— Comparison  of  experimental  and 
theoretical  spectra. 
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L01  ■ 1/2  meter  [87] 

Lq2  “ 100  meters 
z^  ■ 20  meters  [99] 
a ■ 47.5 
8 >0.05 

for  a ratio  of  about  three  orders  of  magnitude  between  the  localized 
and  distributed  turbulence  strengths.  Also  shown  in  this  figure  are 
some  typical  experimental  spectral  data  (Indicated  by  the  points) 
derived  from  the  ALP  program  (run  7 of  flight  72).  The  experimental 
spectra  were  calculated  from  data  obtained  by  directing  a laser  beam 
through  an  aircraft  window  towards  another  aircraft  one  kilometer  dis- 
tant. The  beam  waist,  located  one  meter  Inside  the  window,  had  a 
spotslze  of  20.6  x 10"6  meter.  Flguhe  24  demonstrates  the  ability  of 
the  dirty  beam  model  to  predict  the  low  frequency  breakpoint,  If  not 
the  slopes,  which  are  found  In  the  experimentally  obtained  spectra. 

This  concludes  the  analysis  of  the  temporal  scintillation  spec- 
trum of  a dirty  gausslan  beam. 


D.  Summary  and  Conclusions 

This  chapter  has  consisted  of  two  basic  parts.  The  first  part 
was  devoted  to  developing  a model  of  a dirty  laser  beam,  l.e.,  a laser 
beam  blemished  In  a spatial  sense.  In  the  second  part  the  scintilla- 
tion spectrum  of  a specific  dirty  laser  beam  was  calculated.  The 
analysis  showed  the  scintillation  spectrum  of  the  dirty  laser  beam  to 
be  sensitive  to  turbulence  localized  In  the  vicinity  of  the  transmitter. 
An  example  of  this  localized  turbulence  would  be  the  turbulent  boundary 
layer  near  an  aircraft  fuselage.  The  combination  of  dirty  laser  beam 
and  localized  turbulence  gave  rise  to  the  low  frequency  effects 
described  in  Chapter  III. 

A model  of  the  dirty  laser  beam  was  developed  In  Section  B.  The 
spatial  corruption  of  the  laser  beam  was  modeled  as  arising  from  the 
beam's  being  transmitted  through  a window  whose  transmission  function 
was  unity  except  within  a small  axially  located  region  where  the  phase 
of  the  incident  beam  was  shifted  slightly.  Using  Bablnet's  Principle 
the  equivalent  fields  were  expanded,  within  the  transmitter  plane.  In 
a series  of  eigen-modes  of  the  free  space  wave  operator.  This  expan- 
sion enabled  the  receiver  plane  fields  to  be  written  down  inmedlately 
by  inspection. 
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The  phase  and  diameter  of  the  phase  obstacle  were  assumed  to  have 
been  chosen  a priori  and  an  optimized  expansion  was  developed  which 
minimized  the  integral  square  error  between  the  actual  complementary 
field  and  its  truncated  series  approximation.  This  optimization  scheme 
was  used  to  approximate  the  complementary  field  by  Its  lowest  order 
eigen  mode  so  that,  essentially  the  dirty  beam  consisted  of  the  weighted 
sum  of  two  gausslan  beams. 

Section  C of  this  chapter  dealt  with  a computer  evaluation  of  the 
scintillation  spectrum  of  a particular  formulation  of  the  dirty  gausslan 
beam.  Analyses  were  provided  for  three  different  index  of  refraction 
spatial  spectra:  1)  homogeneous  turbulence  with  a large  outer  scale,  2) 
localized  turbulence  with  a small  outer  scale,  and  3)  a combination  of 
these  two.  The  third  turbulence  distribution  is  appropriate  when  the 
laser  beam  originates  from  an  aircraft.  Experimental  evidence  [87] 
supports  this  assumption  of  a turbulent  boundary  layer  typified  by 
outer  scales  on  the  order  of  one  meter  or  less.  The  results  of  the 
evaluation  showed  that  the  scintillation  spectrum  of  a dirty  beam  dis- 
plays behavior  similar  to  that  of  a clean  laser  beam  with  an  off-axis 
detector.  In  addition,  it  was  demonstrated  that  the  scintillation 
spectrum  of  a dirty  beam  can  display  Increased  sensitivity  to  inhomo- 
geneous turbulence  strengths. 

This  completes  the  analysis  of  the  scintillation  spectrum  of  a 
dirty  gausslan  beam.  Although  the  model  developed  herein  was  a simple 
one,  we  have  shown  it  to  qualitatively  describe  a much  broader  class  of 
situations.  In  practice,  an  actual  laser  beam  may  be  expected  to  con- 
tain a much  more  complicated  spatial  noise  component.  However,  we  have 
demonstrated  that  the  scintillation  spectrum  is  largely  independent  of 
the  actual  structure  of  this  spatial  noise. 

In  the  next  chapter  we  shall  summarize  the  work  presented  in  this 
and  the  preceding  chapters,  and  present  additional  discussion  and  com- 
ments. 


CHAPTER  V 


SUMMARY  AND  DISCUSSION 


A.  Summary  and  Conclusions 

The  effort  described  herein  was  aimed  at  developing  a simplified 
formula  for  temporal  scintillation  spectra  and  applying  It  to  the 
analysis  of  the  scintillation  spectrum  of  a dirty,  l.e.,  spatially 
corrupted,  laser  beam.  The  major  accomplishments  of  this  work  were: 

1)  the  derivation  of  a general  but  very  compact  expression 
for  the  temporal  scintillation  spectrum  of  an  arbitrary 
extended  source, 

2)  Interpretation  of  the  scintillation  spectrum  formula  In 
terms  of  a phase  grating  model  of  the  index  of  refraction 
field, 

3)  asymptotic  and  numerical  calculation  of  spectra  for  a 
gaussian  beam  with  an  off-axis  detector, 

4)  modeling  of  a spatially  blemished  laser  beam,  and 

5)  qualitative  but  not  exact  quantitative  agreement  between 
theoretically  predicted  and  experimentally  obtained 
scintillation  spectra. 

Use  was  made  of  the  so-called  Extended  Huygens  Fresnel  integral 
to  derive  a formula  for  the  temporal  scintillation  spectrum  of  an 
arbitrary  source.  This  expression  was  in  terms  of  the  second  order 
spatio-temporal  statistical  moments  of  spherical  waves.  Expressions 
for  these  statistical  moments  were  derived  by  the  method  of  smooth 
perturbations,  the  weak  turbulence  approximation  was  employed,  and 
various  observations  made,  to  produce  the  final  expression  for  the 
spectrum.  This  final  formula  for  the  scintillation  spectrum  was  in 
terms  of  the  free-space  receiver  plane  fields  of  an  unspecified 
extended  source,  and  in  terms  of  an  extended  receiver.  A phase  grat- 
ing model  of  the  atmospheric  turbulence  field  provided  an  interpre- 
tation of  this  expression. 

The  general  formula  was  used  to  derive  an  expression  for  the  scin- 
tillation spectrum  of  a clean  laser  beam  with  an  off-axis  point  detec- 
tor. A thorough  asymptotic  evaluation  of  this  expression,  valid  for 
arbitrary  beam  Fresnel  number,  was  provided  and  numerical  results 
were  given  from  a computer  evaluation  for  a typical  gaussian  beam. 
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As  a result  of  the  analysis  It  was  shown  that,  by  moving  the 
detector  off-axis  to  the  beam,  an  additional  breakpoint  or  peak 
(depending  upon  the  relationship  of  the  detector  and  beam  axis  with 
respect  to  the  wind  direction)  was  Introduced  In  the  spectrum  at  a 
frequency  below  the  Fresnel  breakpoint. 

The  dirty  laser  beam  was  envisioned  as  arising  from  the  propaga- 
tion of  an  Initially  clean  gausslan  beam  through  a window  and  then 
Into  the  turbulent  atmosphere.  The  window's  transmission  function 
was  assumed  to  be  unity  except  within  a small  axial  region  In  which 
the  Incident  field  was  shifted  slightly  In  phase.  Bablnet's  principle 
was  used  to  calculate  the  equivalent  transmitter  plane  field  which  was 
then  expanded  In  a series  of  eigen  modes  of  the  free  space  wave  opera- 
tor. This  expansion  enabled  the  calculation  of  the  free  space  receiver 
plane  field,  as  required  In  the  general  formula,  by  Inspection  rather 
than  by  performing  the  diffraction  Integral  over  the  transmitter  plane 
field.  The  series  approximation  to  the  transmitter  plane  field  was 
truncated.  A procedure  therefore  was  given  which  optimized  the 
expansion. 

This  modeling  procedure  was  applied  to  the  description  of  a par- 
ticular dirty  beam  and  the  resulting  expression  was  evaluated  numeri- 
cally on  a digital  computer.  The  results  showed  the  dirty  beam  scin- 
tillation spectrum  to  be  sensitive  to  turbulence  In  the  vicinity  of 
the  transmitter  despite  the  fact  that  the  beam  Fresnel  number  was  very 
large.  Specifically,  the  spectrum  displayed  a low  frequency  outer 
scale  dependent  behavior  like  that  of  a clean  gausslan  beam  with  an 
off-axis  detector.  It  was  demonstrated  that  this  particular  model  of 
a dirty  laser  beam  qualitatively  predicts  the  scintillation  spectrum 
of  a much  more  complicated  dirty  gausslan  beam. 


B.  Discussion 

The  general  formula  for  the  temporal  scintillation  spectrum  Is 
of  significant  worth  because  of  Its  simplicity,  flexibility,  and  easy 
Interpretation.  Its  principle  utility  Is  In  predicting  the  spectrum 
for  complicated  extended  source  fields  (for  example  the  field  of  a 
laser  with  cassegralnlan  optics,  or  a spatially  noisy  laser  beam)  and 
for  complicated  receiver  apertures.  With  this  formula  It  Is  possible 
therefore  to  easily  predict  the  scintillation  spectrum  under  realistic 
conditions  rather  than  the  highly  Idealized  situations  which  other 
analyses  presume. 

The  major  limitation  of  our  formula  Is  Its  restriction  to  the 
weak  turbulence  regime.  Thus  one  of  the  suggestions  for  future 
analytic  work  concerns  Its  generalization  to  the  strong  turbulence 
regime.  There  are  two  methods  by  which  this  problem  may  be  attacked. 
First,  expressions  for  the  spherical  wave  statistical  moments  could 
be  found  which  are  valid  within  the  multiple  scattering  domain. 
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Second,  an  heuristic  derivation.  In  the  manner  of  Section  E,  Chapter 
II,  using  a phase  grating  model  of  the  turbulence  field  could  be 
attempted. 


Another  topic  of  further  analytical  effort  could  be  the  consider- 
ation of  the  detector  location  as  being  a random  function  of  time. 

This  Is  equivalent  to  the  laser  pointing  direction  being  a random 
variable.  Under  the  assumption  of  an  appropriate  statistical  behavior 
for  the  detector  location,  It  should  be  possible  to  determine  the 
effect  upon  the  measured  temporal  scintillation  spectrum. 

Throughout  this  paper  we  have  assumed  the  spatial  noise  structure 
of  the  Initial  laser  beam  to  be  non-periodic.  Perhaps  a more  realistic 
model  of  a dirty  laser  beam  should  also  contain  a certain  amount  of 
periodic  spatial  noise.  This  situation  would  be  encountered  when  the 
laser  beam  Is  passed  through  a window  and  Is  multiply  reflected  at  the 
two  surfaces. 

For  future  experimental  effort  we  have  a nuirtber  of  suggestions. 

The  easiest  experiment  to  perform  would  be  aimed  at  verifying 
the  predicted  low  frequency  behavior  of  the  spectrum  due  to  an  off- 
axis  detector.  This  should  be  a simple  method  of  determining  the 
outer  scale  of  turbulence. 

In  case  the  laser  beam  Is  being  aimed  by  a pointer- tracker.  It 
would  be  Interesting  to  calculate  the  temporal  power  spectra  of  the 
servo  system  error  signals. 

Finally,  we  suggest  a closely  controlled  experiment  In  which  the 
source  field  Is  perturbed  In  a known  manner.  The  resultant  data 
should  provide  Information  on  the  routinely  encountered  problem  of 
a dirty  laser  beam. 
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APPENDIX  A 

THE  EXTENDED  HUYGENS-FRESNEL  INTEGRAL 


This  appendix  Is  devoted  to  a discussion  of  the  application 
of  the  Green's  function  technique  to  the  problem  of  propagation  with- 
in a randomly  Inhomogeneous  medium.  In  the  literature  this  general- 
ized Green's  function  technique  Is  commonly  referred  to  as  the 
Extended  Huygens-Fresnel  Principle. 

It  can  be  demonstrated  [100]  that  at  optical  frequencies 
(wavelengths  much  shorter  than  the  Inner  scale  of  turbulence) 
polarization  effects  are  negligible.  We  are  therefore  justified 
In  employing  the  scalar  form  of  the  wave  equation. 


[v2+k2n2(*)]E(K)  ■ 0,  (A-l ) 


where  the  Index  of  refraction  Is  a random  function  of  position.  We 
now  propose  the  existence  of  a Green's  function  such  that 

[72+k2n2(lT)]G(lU?')  ■ -4ir«(|ff-lf'|).  (A-2) 


By  left-multiplying  Equation  (A-l)  by  G(O')  and  Equation  (A-2)  by 
EOT),  subtracting  the  two  resulting  equations,  and  Integrating  over 
the  scattering  volume  we  obtain 


d'R[G(U'  )v2E(K)-E(*)v2G(ff,K' )] 

V 


■ 4w  d*  E(ff)«(|ff-ff'|)  - 4irE(JT' ) . (A-3) 

JV 


Figure  1 Illustrates  the  geometry  of  Interest.  Point  P0  denotes  the 
source  point  and  Pi  the  field  point,  £ and  I*  the  surfaces  bounding 
the  scattering  volume  v and  the  outward  directed  unit  normal  to 
the  surface. 

If  the  first  and  second  partial  derivatives  of  G and  E are 
single  valued  and  continuous  on  £ and  £'  we  can  make  use  of  Green's 
second  Identity,  which  when  applied  to  Equation  (A-3)  yields 
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E(*'>  - W/6H-  Ets)ds- 


(A-7) 


where  in  the  interest  of  notational  simplicity  we  have  dropped  the 
explicit  dependence  of  E and  G upon  and  IT. 

For  the  field  of  a point  source 


E(R)  = S. 


(A-8) 


the  partial  derivative  normal  to  the  surface  z is 


3E  _ rAv 

IN  “ H k rJ  IN  • 


(A-9) 


If  we  assume  R » 1/k,  i.e.,  the  surface  z is  many  wavelengths 
removed  from  the  point  source  at  pQ,  then  we  have 


If  *,kE<R>  If  • 

Using  the  notation 


3R  - _D  . z 
aN  R z » 


(A-10) 


(A-ll) 


in  which  the  circumflex  accent  denotes  a unit  vector,  we  have 


% -ikE(R)R  • z . 


(A-12) 


We  now  claim  [41]  that  the  generalized  Green's  function  can  be 
expressed  as 


G(R.R') 


,ik|I-R'  N(R.R') 


(A-13) 


where  * is  a complex  function  arising  from  propagation  through  the 
inhomogeneous  medium.  The  perturbation,  ip,  to  the  field  of  a 
spherical  wave  is  obtained  by  solving  the  scalar  wave  equation 
(Equation  ( A-l ) ) for  the  Green's  function,  G,  rather  than  for  the 
field,  E.  However,  for  the  purposes  of  this  discussion,  we  need 
only  to  bound  the  magnitude  of  i|j.  This  task  is  easily  accomplished 
by  inspecting  \p  within  the  geometrical  optics  domain  (/xT  < iD). 
Under  this  condition  we  find  that  [41] 


♦(*.*•)  3 ik  | (s)ds 


(A-14) 


where 


s = |R-R'  | 

and  we  have  used  Reynold's  convention  [102]  on  the  Index  field 
n(R)  3 n-jO?)  +<  n(R)>  3 n^R)  + 1. 

The  generalized  Green's  function  as  given  by  Equations  (A-13)  and 
(A-14)  can  be  shown  to  satisfy  the  homogeneous  wave  equation  (A-l) 
if  the  term  proportional  to  vn]  is  neglected.  (The  polarization 
term  [100]  in  the  wave  equation,  which  is  a function  of  vni , was 
neglected  in  obtaining  Equation  (A-l).) 


Now  the  partial  derivative  of  the  Green's  function  is  given  by 

(A-l 5 ) 


«■•[(« -*)*♦«] 


For  s » 1/k  Equation  (A-l 5)  is 

aN*  1kG|w  + nr affj  * 

Within  the  geometrical  optics  domain 


1 

TiT 


Since 


tlR  = j n-j  (s)ds  - |n  |s  { n!  (s)ds 

1 3*  3S  / » 

T73N  3N  nl(s)- 

17  ~ini (s)l  K< 


3 S 

an<*  3N  ^ ^ * Equation  (A-l 5)  becomes 
1^  * -1kG(R,R' ) s • z, 


(A-l 6 ) 


(A-l 7 ) 


(A-l 8 ) 


(A-l 9 ) 
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where  we  have  defined 


r' 


(A-12)  and  (A-19)  the  expression  for  the  field 

G(R,R')E(R)(R-z  - s-z)dA.  (A-20) 


3s  _ : ; 

8N  ' _s 


With  Equations 
at  becomes 


E(R') 


If  we  make  the  paraxial  approximation,  i.e.,  restrict  our  attention 
to  small  angles  with  respect  to  the  axis  of  propagation,  we  can 
approximate  the  obliquity  factor  [103]  and  the  Green's  functions 
respectively  by 

A A A A 

R*z  - s*z  ^ 2 


and 


G(R,R' ) % 


ik| R-R*  |+if>(R,R' ) 


With  these  approximations  Equation  (A-20)  becomes 


E(R')  = E(R)elklR_R'  l+^(R»R‘ ) dR  . (A-21 ) 

Equation  (A-21)  is  the  final  result  of  this  appendix. 

Although  in  the  preceding  discussion  we  restricted  our  attention 
to  the  case  in  which  ^ was  purely  imaginary,  it  is  in  general  complex. 
This  approximation  was  merely  a device  to  enable  the  neglect  of  the 
term  Involving  the  partial  derivative  of  i/>.  The  actual  expression 
for  if>  is  obtained  by  solving  the  scalar  wave  equation  for  the  gener- 
alized Green's  function.  This  Green's  function  is  simply  the  field 
(which  has  been  perturbed  in  phase  and  magnitude  due  to  propagation 
through  the  turbulent  atmosphere)  Of  a point  source.  An  approximate 
expression  for  the  complex  phase  perturbation,  ip , is  derived  by 
Tatarski  [52]. 


Equation  (A-21)  is  actually  more  general  than  it  would  seem 
from  our  derivation.  This  formula  states  that  the  field  at  R'  is 
given  merely  by  the  convolution  of  the  source  field,  at  the  surface 
R',  with  the  generalized  Green's  function.  Stated  another  way,  the 
field  at  R'  is  due  to  the  linear  superposition  of  the  elementary 
wavelets  originating  at  the  plane  R. 


APPENDIX  B 


DERIVATION  OF  THE  SPHERICAL  WAVE  FIRST  ORDER 
PERTURBATION  SPATIAL  SPECTRUM 


In  this  appendix  we  shall  supply  the  missing  steps  between  Equations 
(45)  and  (46)  of  Chapter  II.  The  beginning  point  is  Equation  (45)  which 
relates  the  spatial  spectrum  of  the  complex  phase  perturbation  to  the 
spatial  spectrum  of  the  index  fluctuations  v; 


F 


(■ejxj+ic£yj) 


d^Uj.iCg.Ljt)  = 


4 fr-efty 


iUjx+K^y) 


'(xi-x)2+(yj- 

8ik[  ?nr-7) 


:j-x)2+(y]-y)2  (x-xi)2+(y-yi)2  (xj-xi)2+(yj-yi)2 

2TlT)  + 2z  ■ K 


(B-l ) 


If  Eq.  (B-l)  is  now  multiplied  by 
dK]dtc2  xi  +<2yl  ^ 


(B-2) 


and  integrated  over  x,  xj , y,  and  y-j  from  -»  to  °°,  the  left  hand  side  is 


L.H.S.  = | dxif  dy-j  |j  d«(<),<',L;t) 


i (xjx^+K^y^ )-i (icjxj+Kjyj ) 


Interchanging  orders  of  integration  gives 


(B-3) 


00  C • I 

L.H.S.  = dK1d<2  ||  d*(<j,K^,L;t)ii^  | dxj  e™1'^  J 

—oo 

f , r , (k2-ic2 

* Ujjyi  e ’ J 


(B-4) 


which  is 


w 

L.H.S.  = dK1dK2  ||  d$(K^,K^,L;t)6(Kj-K1)6(K^-tc2)  (B-5) 


or  finally 


L.H.S.  = d$(K-j  ,<2,L;t) 


Equation  (B-l)  is  therefore 


(B-6) 


2 L 

d$(Kl,K2,L;t)  = d^d^  | ^(0  ||  dv(<j  ,K2,z;t) 


x [ fdxdy 

4/ 


( i ( K.x+K  _y ) -i  (k-i  X-.  ttcpy-j ) 

dx'dy'  e 1 ^ 1121 


1 1 


-OO  _00 


ik 


x e 


(xj-x)2+(yj -y)2  (x-x- )2+(y-y, )2 

+ 2z 


(xj-xi)2+(yj-yi): 


2L 


1 


(B-7) 


The  x,  x\  integrals  only  are 
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ruj-x)2 

i r r Kieix-K-x.')  Ik  1 2'( L - zT 

h f * Ld*i  e 1 1 • 


(x-X^' 
+ — 2^— 


(xi- 

“TU 


h j>  f dxi 


i ( tc-j'  x-<-j  x-J ) 


(B-8a) 


') 

ik  ■ 
x e 


x.j2-2xjx+x2  x2-2xxi+x2  x^2-2x^xl-+xi2  "| 
2(L-z)  + Tz  2L  I 


(B-8b) 


r xj2  xl2~2xixj] 

-1  ic] xi  Ik  2(L-Zy  ' 2L 
dxj  e e L J 


(B-8c) 


'»  iic-Jx  ik 
x dx  e e 


O p 

-2x-jx+x  x -2xx. 
2(L-zT  ' + 2z 


The  exponent  of  the  x integrand  only  is 


.2  r 1 k . ik  1 o„ 

[s + 2* 


ikx^  iK-j  ikxi 


2[L=I)  ‘ T"  + 2z 


(B-9) 


This  Is  of  the  form  Ax2-2xB.  Completing  the  square  and  performing  the 
x integration  gives  the  following: 


[*  dx  e*x2-2»B  - e-62Mf  dx  eA<*-B'A>2 
' -00  -00 


(B-10) 


The  x,  Xj  Integrals  are  then 
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2 v 2 


r xi2  xi2-2xixi“ 

-ikz(L-z) 

r v 1 '-1 

xi  * _?1  “ll 

[2(L-z)  2L 

e 2L 

z L-z  kj 

(B-lla) 


,2 


^TkL 


i (L-z)x.k ' iz(L-z)k1‘  ,/  , z\ 

1 r , •1xl(KT'cl  l) 

2tt  J dxi  e V 
4 —00 

(B-llb) 


H [ r 2kL 


,2 


i (L-z)xiK-j  iz(L-z)Kj‘ 


It  is  easily  seen  that  the  y,  y-j  integrations  give  similar  results; 

1 (L-zJy^Kg  i z(L-z) <2^ 

k e—  6("2  - z "2)  (W» 


so  that  Equation  (B-7)  is 

2 l 00 

d<*>( <i  ,<2»L;t)  = dic-jd^  2tT  | 

^ “00 

X (if  P^1]  5(‘i  • I ‘l)  S(‘2  - 2 x2) 

- 1Z2k[Z>  (t]2+t22) 


x e 


(B-13) 


Performing  the  *'  integrations  then  gives  the  final  answer; 


d*(Kl,K2,L;t)  = ik 


x e 


iXkii-.-  k 

z ri  2k  z 


(B-14) 


which  is  simply  Equation  (46)  of  Chapter  II. 
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APPENDIX  C 

THE  WEAK  TURBULENCE  APPROXIMATION 


This  appendix  Is  devoted  to  a description  of  the  conditions  under 
which  we  are  justified  In  making  the  weak  turbulence  approximation  as 
It  Is  employed  In  Section  D of  Chapter  II. 


The  weak  turbulence  approximation  requires  that 


M.IbI  « 1 


(c-i) 


where 


° “ T [Dw(^rr2'°i°)  + Dw(r3_^4»0i°)3 


6 * - 7 tvyy'i^'1 + 

- - Vy^4iri^21T)] 

+ 2 [Cx(r1-r3;rj-r^;T)  + Cx(r2-r4;rj-r^;T) 

+ 1CxS<V?3;7i-f2iT)  " 


(C-2) 


xSv  2”  4’  1 2’ 


(C-3) 


DISCUSSION  OF  a 


From  the  results  of  Chapter  II,  Section  C,  viz.  Equation  (76), 
we  can  express  a as 


a = 2k* 


r.  f.1/2  f2z  ,1  ,2(l-z) 

d*  i dz  dp  + dz  dp 

JJ  |Jo  Jo  J 1 /P  Jo 


1/2  J0 


i (1-z)k- (r1-r2)  1(1-z)k* (r3-r,) 

x Fn(*c-| ,k2»pl)  1 - e + 1 - e 


•S*  A 


Since  the  p integrand  Is  non-negative  (Fn  is  a power  spectrum)  we  have 


a < 2k2L2 


fl"7 


i(l-z)f(rrr2)  1(l-z)K-(r3-r4) 

1 - e 1 6 + 1 - e J * 


x Ldp  Fn(VK2*pL) 


(C-5) 


Using  the  relationship  between  the  two  and  three  dimensional  index  of 
refraction  spectra. 


j dp  Fn(K]  *<2»pL)  = 4>n(Kl  »*c2»°) 
0 


and  making  the  change  of  variable 


u = 1 - z 


(C-6) 


(C-7) 


gives  for  Eq.  (C-5) 


°o  ^ 

2irk2L  ||  die  ♦p  ( <1  **c2  »0)  | du 


luic*  (r-|-r2)  1u<*(r3~r4) 


+ 1 - e 


(C-8) 


Changing  the  k variable  of  integration  to  polar  coordinates,  assuming 
isotropic  turbulence,  and  performing  the  angular  integration  yields 


ot  < 4ir^k2L  | die  k #n(<)  | < 


X [1  - JQ(uiep^2)  + 1 " JQ(uiep34)] 


(C-9) 


where 
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and 


p12  * |rr  r2| 


P34  * I r3"  r4l 


Now  for  the  von  K£rm£n  form  of  the  index  of  refraction  spectrum. 


♦n(0  * -033  C2[k2+(1.077/Lo)2]'11/6 


(C-10) 


we  have 


j^d<  *n(K)[l-J0(uKp12)]  * 

.033  cl  jj^dic  k[k2+(1.077/Lo)2]'11/6 

d<  k[k2+(1.077/Lo)2]‘11/6  J0(u<p12) 

| 

The  right  hand  side  of  this  equation  is  easily  put  into  the  form 
.033  C*  (Jtf)'5/3  {f  «x  *(Al) -,1/6 


(c-11) 


y"1^2  j dx  x^x2*!)-11'6  JQ(xy)  v'xy 


(C-12?) 


where 


(C-12b) 


up-i  p 

y = 1.077 

0 

and  we  have  made  use  of  the  change  of  variable  x = kLo/1.077. 

The  first  integral  in  the  above  expression  yields  a Beta  function  [73] 
and  the  second  integral  is  a Hankel  transform  [104].  Performing  these 
integrations  results  in 
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m 

|0dle  K ®n^'c^1"Jo^u'cp12^  * 
•033  cj;  (W6/3  * 


5/6 

, . * 5/6 

1 - 2 575~ 

2 r (5/6) 


(C-13) 


where  Kc#g  Is  the  modified  Bessel  function  of  the  second  kind  and  order 
5/6.  We  require  the  separation  to  be  much  smaller  than  the  outer  scale, 

p12  <<:  Lo  * 

so  that  the  right  hand  side  of  Equation  (C-13)  can  be  approximated  [105] 
as 


.033  c; 


!W 


■5/3 


r(5/6)  r(l/6) 
2r(ll/6)  r(ll/6) 


r 


or 


5/3 

•033  C2 


With  these  results,  inequality  (C-9)  gives  for  a. 


a < 4,2k2L(.033)C2 


w 


5/3 


(^5/3- 


du  u 


5/3 


or  finally 


“ < 10 


*033,C2  .2k2L  » [(^)  * (Sif  ] 

The  separations,  p^j,  are  bounded  by 


(C-14) 


(C-15) 


(C-16) 


, < XL 

*1J  * IT 
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where  iQ  is  the  inner  scale  of  turbulence.  For  a Helium-Neon  laser 
(x  * .6328  x 10"6  meter),  a one  centimeter  inner  scale,  and  a one 
kilometer  range  we  obtain 

a < 5.384  x Id'6  * p^3)  4 lo'5  cj;  . 


(C-17) 


Therefore  for  1015  C2  < 1 we  can  employ  the  weak  turbulence  approxima- 
tion. 

DISCUSSION  OF  e 

With  Equations  (72),  (74)  and  (76)  in  the  text,  6 Is  expressed  as 

, “_f.l/2  ,2z  .1  ,2(1 -z)  ] 

6 = 2kzl/  j jdtc  -«  dz  j dp  + j dz  J dp  |>  Fn(*i  ,lc2»pL' 


x e 


10  J0 
itc*[z(r|-r^)+Tv] 


‘ 1 (1-z)k- (r1-r4)  i(l-z)K*(r2-r3) 

e + e 


i(l-z)»c*(r-| -r3)-i2y  i(l-z)ic*  (r2~r4)+i2y 
e - e 


(C-18a) 


where 


y = 


(z-z2-p2/4)  st- 


(C-18b) 


Now  e is  maximum  for  zero  time  lag(T=0)  and  zero  receiver  plane 
separation  (rj-r2=0)  so  that 


? ? f7  _ ff1/2  f2z  f1  f2(1_z)  1 , 

<_  2kzl/  J J die  i J dz  J dp  + ^dz  dp  - Fn(iO|  .K2,pL) 


0 '0 
1(l-z)lT-r14  i(l-z)7*r23 


+ e 


i(l-z)K*r-j3-i2Y  i(l-z)K*r24+12Y 


- e 


- e 


24' 


(C-19) 
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r 

j. 


where  r^j  * r1-rj* 

Further,  B Is  maximized  by  setting  the  separations  Fjj  equal  to  zero: 


B $ 4kzL 


x < 


1 - cos 


,7-  ffV2  f2z  fl  f2(l-z)  1 
2,2  die  y dz  dp  + dz  dpf  F (k,,ic2,pL) 

ll  J0  J0  J 1/2  J0  I 
[(z-i2-cz/4>  y] 


Little  error  Is  incurred  by  dropping  the  p dependence  of  the  cos 
function  and  extending  the  upper  limits  on  the  p integrals  to  infinity 
(see  Appendix  D)  to  obtain 

P £ 4nk2L  JJd<  j^dz  *n(i<i  **2.0)  -f 1 - cosp^~z)L  ic2ll  . (C-21) 


Changing  the  k variable  to  polar  coordinates,  assuming  isotropic  turbu- 
lence, and  performing  the  angular  integration  gives 


e £ 


.1 

8ir2k2L  die  k 4 (k)  f dz 

Jo  " Jo 


1 - cos 


gg-*?L  k2 

k 


(C-22) 


With  the  identity 


1 - cosx  = x 


fl 

du  sin 

* n 


UX 


(C-23) 


inequality  (C-22)  becomes 


B $ 8ir2k2L  j^dx  k *n(ic)  j^dz  <2 j 

!]  • 


( du  sin  u K2 

Jq  L k 


(C-24) 


or 


133 


0 % 8iT2k2L 


du 


i , x . r z0-z)l  21 

die  k 4>n(ic)sin  N k k J 


(C-25) 


For  the  von  KarmSn  index  of  refraction  spectrum,  as  in  Equation  (C-10), 
we  find 


* * -033  Cp  87t2k2L  j^dz  ["(Lk^kj  |Qdu 
x | d<  k3[k2+(1.077/Lo)]'11/6  sin  [u  z^'z^-  <2J  . (C-26) 


A stronger  inequality  results  if  the  outer  scale  is  set  to  infinity 
(1.077/Lo=0); 


e < .033  i eA2i 


I'd,  I du 


■]fr 


f die 
J n 


<-2/3  sin 


zlhzjL  K2l 


(C-27) 


or 


0 < .033  C2  4ir2k2L 


rl 


f'dz  [^1  fdu 

Jn  L K J J c 


x dx  x-5/6  sin 
J0 


0 >-  " -1  J0 
z(1-z)L 


r ziLiiL  1 

U k X 


(C-28) 


where  we  have  made  use  of  the  change  of  variable 
2 

x = KC  . 

Performing  the  x integration  [106]  yields 
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6 < .033  Cp  4»‘k*L  j^dz  |£Llj3ikj  JQdu 


§■  ’I 

[u 

r0  - 1 

COS  ^ 

1 . 5 

2 6 v 

> 

(C-29) 

or 

6 < .033  cn2  4ir2k7/6L1 1/6  r(l/6)  cos  (fj) 

X f’dz  z^O-z)5'6  f'du  u’1/6  . (C-30) 

Jo  Jo 

The  u integral  gives  the  value  6/5  and  the  z integral  is  a Beta  function 
[73],  Finally  we  have 


3 < .022  C2  w2k7/6L11/6 


9 

5 


cos 


5tt  r(l/6)r2(ll/6) 

12  r(8/3) 


(C-31 ) 


From  Tatarski  [21],  the  log-amplitude  variance  for  a spherical  wave  with 
an  infinite  outer  scale  of  turbulence  is 


a2  = .124  C2  k7/6L11/6  (C-32) 

A * * 

so  that 


0 < 4a2  . (C-33) 

Inequalities  (C-17)  and  (C-33)  constitute  the  results  of  this 
appendix.  These  expressions  state  that  we  are  justified  in  employing 
the  weak  turbulence  approximation  in  the  development  of  Chapter  II  if 

5.384  x 1013  C2L  (p125/3  + p345/3)  « 1 
and 


4p?  « 1 


APPENDIX  D 

DEPENDENCE  UPON  LARGE  OUTER  SCALES 


In  this  appendix  we  discuss  the  approximations  Involved  In  going 
from  Equation  (88)  to  Equation  (89)  In  Chapter  II. 


The  scintillation  spectrum  Is  expressed  by  Eq.  (88). 

r2(l-«) 

lo  jo  n/2  Jo 


, o rl/2  f2z  fl  r2(l-z) 
SD(u>)  ■ 2kzLz  J dz  dp  + dz  dp 

p I Jo  JO  *1/2  Jo 

00 

x JJdK  Fn(<]  ,K2 >pL)  [2ir6(u-ic*17)] 

"00 

x |H(z,p,7,w)  - H*(z,p,-7,w) |2 


(D-la) 


where 


H(z,p,k,w)  ■ 

-1  (1-z-pV  4)Lk2  __  , 

e « Jdrw(F)  •1l',r  E (f  - E*(r) 


(D-lb) 


We  will  show  that  under  a reasonable  set  of  circumstances,  the 
p dependence  of  the  H function  can  be  suppressed  and  the  limits  on 
the  p Integrations  extended  to  Infinity.  First  we  shall  demonstrate 
the  approximations  for  high  spatial  frequencies  and  then  for  low 
spatial  frequencies. 

To  simplify  the  discussion  which  follows,  we  adopt  the  following 
notation: 


H+  = H( z,p,k,w) 


= H(z, p,-k,w) 


(D-2a) 


(D-2b) 


The  square  modulus  term  In  Equation  (D-la)  Is  then  given  by 
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|h+-h;| 


dFw(F)eiK’r  E^F  - E*(r) 


drw(r)eiic‘r  E(r)E*^r  + 


- e 


-i  (l-z-p2/4)^- 


2 r 


drw(F)  eiK'r  E^r  - k)  E*(r) 


dFw(r)e"lK’r  E*(r)E^r  + <) 


- e 


i(l-z-p2//4) 


2 . 


drw(r)e~iK’r  E*^r  - E(r) 


x 


drw(r)e’K  r E(r)E*^r  + ^ 


(D-3) 


For  high  spatial  frequencies  (k  » 2ir//xT)  the  exponential  terms 
involving  p oscillate  so  rapidly  with  respect  to  the  other  terms  in 
the  expression  that  we  may  replace  them  by  their  average  integral 
value,  zero.  Within  the  high  frequency  region  then 

|H+-H?|2  = |H+|2  + |H_ 1 2 , (D-4) 


and  the  dependence  upon  p disappears  from  this  term.  The  onset  of 
this  region  is  at  roughly 

L<2  „ ? 

T-  ^ 2* 

that  is 


K ^ 27T//XL 


(D-5) 


However,  under  the  assumption  of  isotropic  turbulence,  the  two  dimen- 
sional index  of  refraction  spectrum  has  a breakpoint[57]  at 


k-pL  'v*  2tt 


(D-6) 


beyond  which  the  spectrum  decreases  rapidly.  So,  for  high  spatial 
frequencies,  Fn  is  very  small  beyond  the  point 


2tt_ 

/ILP 


'V  2tt 


or 


p 


«<  1 


(D-7) 


Therefore  we  may,  with  little  error,  extend  the  upper  limit  on  the  p 
integrals  to  infinity  to  obtain 


r1/2  f 

2z 

fl 

r2(l-z) 

< 

dz 

dp  + 

dz 

dp 

0 J 

0 

1/2  J 

0 

Fn^Kl  »tc2’p*-) 


x |H(z,0,k,w)-H*(z,0,-k,w)|2 

= f *n(Kl»K2»°)  |H(z.0,7,w)-H*(z,0,-k,w)|2  , (D-8) 

J0 

where  we  have  made  use  of  the  relationship  between  the  two  and  three 
dimensional  index  of  refraction  spectra; 


— 

L 


(D-9) 


Equation  (D-8)  is  the  desired  result  of  this  appendix.  We  have 
justified  this  approximation  for  high  spatial  frequencies  (k»2t r//xL). 
It  will  now  be  demonstrated  that  Equation  (D-8)  is  valid  at  low  fre- 
quencies as  well. 

Now  recall  that  the  two  dimensional  index  of  refraction  spectrum, 
Fp*  is  small  beyond  the  point  <pL  'v  2ir.  Therefore  the  frequencies  for 
which  Fn  is  significant  are  given  by 


kpl  Z 2n  (D-10) 
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*• 

■ 


or 


or 


(<pL)2  < ( 2tt) 2 


l££i£k  < » ! « * 

8k  ^ 4 I. 


(D-ll) 


The  left-most  term  in  Inequality  (D-ll)  is  merely  the  p dependence  of 
the  H function  (Equation  (D-lb)).  Since  the  p dependent  term  of  the 
H functions  is  much  less  than  -n  within  the  region  where  Fn  is  signifi- 
cant, i.e.,  low  spatial  frequencies,  we  can  drop  this  term  altogether. 
The  H function  therefore  becomes 


H(z,k,w)  = e 


■i (I-z)Lk2 
2k 


drw(r)  e 


:(f  - 


E*(r) 


(D-12) 


For  low  spatial  frequencies  such  that 

k « 2ir//\i 

the  field  translation  term  in  the  expression  for  the  H function  is 
bounded  in  the  following  manner: 

L(1-z)k  L 2tt  r-^~ 

k/ff=Ar  ' (D‘13) 

For  well  behaved  (free  space)  receiver  plane  fields,  i.e.,  those  vary- 
ing slowly  over  a lateral  distance  on  the  order  of  a fraction  of  the 

term  in  Equation  (D-la)  can  be 

( D-l 4a ) 

, (D-14b) 


Fresnel  zone,  /XL,  the  square  modulus 
wri tten 


|H+-H*|^  % 4C2  sin2 


(1-z) 


u£l 

2k  J 


or  further 


* 4C2  (1-Z)  ^ 


where 
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A V*  ftl  Ul/  «*' 


w(r)  | E(r)  | ‘ 


(D-15) 


Under  most  cases  of  interest  (receiver  pupil,  w,  smaller  than  a 
Fresnel  zone,  /XL)  the  C function  is  a very  weak  function  of  < for 
k<2tt//xT.  Hence  the  low  frequency  dependence  of  the  modulus  squared 
term  is  k*.  This  fourth  power  dependence  upon  spatial  frequency 
strongly  suppresses  any  contribution  to  the  scintillation  spectrum 
due  to  low  spatial  frequency  atmospheric  inhomogenieties.  For  the 
Kolmogorov  (infinite  outer  scale)  model  of  the  index  of  refraction 
field  the  low  frequency  dependence  of  F is  k"8/3.  Even  for  this 
model  the  *4  dependence  of  the  square  modulus  term  will  cause  the 
scintillation  spectrum  to  be  highly  insensitive  to  low  spatial  fre- 
quency index  of  refraction  Inhomogenieties. 

With  the  foregoing  arguments  the  scintillation  spectrum  is  expressed 


Sp(o>)  = 4ir2k2l  JjdK  ,k2»0)6(u)-k*v) 

-oo 

x f dz  |H(z,7,w)-H*(z,-ir,w) |2  , 

J n 


(D-16a) 


where 


H(z,k,w)  = 


APPENDIX  E 


DIFFRACTION  OF  AN  ARBITRARY  FIELD  BY  A PHASE  GRATING 


Within  this  appendix  we  develop  a mathematical  description  of  the 
fields  diffracted  by  a weak  phase  grating. 

We  assume  the  same  geometry  as  In  Fl£ure  8 of  Chapter  Ilia  phase 
grating  of  period  2*/ I* I and  orientation  k/|7|  located  In  some  plane 
z such  that  0<z<L.  The  weak  phase  grating  Is  defined  In  terms  of  Its 
transmission  function 


1ecos(<"  • 7 ) 

T(<  ,f)  ■ e 


(E-l) 


where  |e|«l. 

Keeping  only  the  first  two  terms  In  the  series  expansion  of  Equation 
(E-l)  yields 

_ _ _ _ e IF  *r  e -1iT  *7 

T(k  ,r)  % 1 + 1e  cos(k  *r)  ■ 1 + 1 2 e +1  2 e (E-2) 


Now  assume  an  arbitrary  field  to  be  Incident  upon  the  grating  from 
the  left.  At  z-z"  the  field  Is 

E(r)|^_  , (E-3) 

and  at  z=z+  the  field  Is 

1 Inir  *7 

£(r)  + a l , cn  e (E-4) 

z n--l  z- 


where 

In  I 

Cn  = (1  e/2)  1 (E-5) 


The  receiver  plane  field  Is  determined  by  performing  the  Fresnel 
diffraction  Integral  over  the  field  within  the  plane  z=z*; 
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E(r' 


') 

L 


jML-zir  ^ 


ik 


L 2tt(L-z)  j 


| dr  E(r) 


2(  L-z) 


(r2-2r*r'+r'2) 


(E-6) 


Combining  Equations  (E-4)  and  (E-6)  yields 


E(r') 


. Jk(L-2)  f -ik  1 i 


L 

ik 


x e 


2(L-z) 


L2tt(L-z) 

!r2-  2r  •[?  - nikil  -j 


ikr'2 

-1  e2(L'z)  l Cn  [dr  E(r) 
J n=-l  J 


(E-7) 


Symbolically  performing  the  integration  gives  the  receiver  plane 
fields 

. 2 


ik(L-z)  n i n<  *r'  ~in  fL~z)  2 / nfL-zl  \ 

1 e ^ Y’-’k'  7)-  (E-8) 


n=-l 


This  expression  for  the  principal  diffracted  fields  (minus  one,  zeroth, 
and  plus  one)  of  a weak  phase  grating  is  the  result  of  this  appendix. 
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APPENDIX  F 

DETAILS  OF  THE  DIRTY  BEAM  MODEL 


The  objectives  of  this  appendix  are  twofold.  First  we  shall 
demonstrate  that,  with  little  loss  of  generality,  we  can  assume  the 
phase  obstacle  to  be  axially  located  within  the  beam.  Then  we  pre- 
sent an  outline  of  the  optimization  of  the  eigen-mode  expansion  of 
the  complementary  field. 

We  restrict  our  attention  to  the  situation  in  which  the  receiver 
is  in  the  far-field  of  the  main  laser  beam,  (Fresnel  number  much 
larger  than  unity)  and  the  detector  is  near  the  beam  axis  (r0«w(L)). 
Under  these  conditions  the  scintillation  spectrum  of  the  main  beam 
alone  is  that  of  a spherical  wave  and  is  therefore  insensitive  to 
the  detector  location.  Now  consider  the  case  in  which  the  receiver 
plane  field  consists  of  the  field  of  the  main  beam  plus  the  diffrac- 
tion pattern  of  the  phase  obstacle.  As  long  as  the  detector  is  near 
the  axis  of  the  main  beam  the  scintillation  spectrum  (of  the  dirty 
beam)  will  be  dependent  only  upon  the  relative  positions  of  the 
detector  and  the  phase  spot  diffraction  pattern.  Therefore  we  can 
assume,  without  loss  of  generality,  that  the  phase  obstacle  (and 
hence  its  diffraction  pattern)  are  axially  located  with  respect  to 
the  main  beam. 

Now  in  performing  the  actual  eigen  mode  expansion  of  the  com- 
plementary field,  there  are  two  points  with  which  we  must  concern 
ourselves:  1)  The  phase  curvatures  (within  the  z=0  plane)  of  the 
main  and  complementary  fields  must  be  identical  (so  that  the  mode 
expansion  will  indeed  converge  to  the  actual  complementary  field). 

2)  Since  we  wish  to  approximate  the  complementary  field  by  only  a 
finite  number  of  terms,  we  must  insure  a close  approximation  to  the 
actual  field,  i.e.,  optimize  the  expansion. 

The  phase  curvature  of  the  complementary  field  is  contained  in 
the  weighting  function  of  the  eigen  mode  expansion,  (see  Equations 
(200)) 


kM/2L  , 2a  2. 
l+iM(z+z-|  )/L  (x  y } 


(F-l ) 


Comparison  of  this  weighting  function  with  that  of  Equation  (194) 
shows  that  for  equal  phase  curvatures  (at  the  plane  z=0)  we  must  have 
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Im{a(z0)>  ■ lm{0(zi ) } 


that  Is 


g0(N/L)2  _ Z] (M/L)2 

l+(z0N/L)2  l+^M/L)2 


(F-2) 


Our  optimization  procedure  Is  based  on  the  minimization  of  the 
Integral  square  error  between  the  complementary  field  and  Its  truncated 
series  approximation.  The  error  Is  defined  as 


OB 

eL  5 II  |Ec(x.y,0)-E[(x,y.0)|2  dxdy 


(F-3) 


where 


C L L 

EL(x.y.2)  - I JE  %„  ^(x.2+z,  )a„(y. Itii 


m«0  n»0 


(F-4) 


By  making  use  of  the  orthonormal Ity  of  the  mode  functions  and  Equation 
(198)  for  the  complementary  field  we  find  that  Equation  (F-3)  can  be 
expressed  as 


w 

b L * dxdy  I D^(x-xs.y-ys)E(x,y,0)| ' 


L L 


’ I-  l-  lYmnl 


m»0  n*0 


(F-5) 


The  error  Is  then  minimized  by  taking  the  partial  derivative  of  el  with 
respect  to  the  complementary  beam  Fresnel  number: 


3 

iff 


L L o 

I l l*mJ  = 0 

m*0  n*0 


(F-6) 
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The  entire  problem  now  reduces  to  that  of  solving  two  simultan- 
eous equations  (Equations  (F-2)  and  (F-6) ) for  the  complementary  beam 
waist  spot-size  Wj  and  locations  z^. 

We  now  proceed  with  an  example  of  the  use  of  the  foregoing  com- 
plementary field  expansion.  In  the  interest  of  mathematical  simpli- 
city the  complementary  field  Is  approximated  by  only  Its  lowest  order 
eigenmode.  The  expansion  coefficient  of  this  mode  Is  given  by 

oo  . 

y00  = {}  dxdy  EC(x’y’°)  9S(x,zi)g*(y,z1)  . (F-7) 


As  has  been  justified  previously,  without  loss  If  generality  we 
can  assume  the  phase  object,  described  by  Dfi,  to  be  axially  located. 
That  is 


Da(x.y)  = 


A eiu  ; J x2+y2  <_  e/2 
1 ; elsewhere 


(F-8) 


Employing  Equations  (194),  (198),  (200),  (F-7),  and  (F-8)  we  obtain 


y 


oo 


WqW^I-A  eiy] 


a(zQ)B*(z1) 

a(z0)+B*(z1) 


f -[a(z  )+B*(z1)](6/2)2 
x <1  - e 0 1 


(F-9) 


Application  of  Equation  (F-6)  to  y00  yields,  after  some  lengthy 
algebraic  manipulations 


x = 1 ^ [l-2xy(T§)]  - 1 

where  we  have  defined 

x . ,r<11) 

“r(2o) 


(F-lOe) 


(F-lOk) 
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y = ctr(z0) 


(F-lOc) 


In  Equations  (F-10)  x is  seen  to  be  the  square  of  the  ratio  of  the  main 
beam  to  complementary  beam  spotsizes  (at  z=0)  and  y the  square  of  the 
ratio  of  the  spot  radius  to  the  main  beam  spotsize.  Solution  of 
Equation  (F-lOa),  which  at  first  glance  appears  a formidable  task,  is 
easily  accomplished  by  an  iterative  technique; 


xn+l  = y An  l"2yx 


in  which  we  let 


■(S3  ■ 


(F-lla) 


= u-1 


(F-llb) 


Equation  (F-llb)  for  the  initial  guess  for  x states  simply  that  (at 
z=0)  the  complementary  beam  spotsize  will  be  on  the  order  of  the 
radius  of  the  phase  obstacle. 

From  Equation  (D-lOb)  we  see  that 


Br(z] ) = X a (z0) 


(F-12) 


and  from  Equation  (D-2) 


®i  (zi ) - 


(F-13) 


Since  x is  known  (from  the  solution  of  Equation  (F-lOa) ) and  a(z0)  is 
known,  we  have  two  simultaneous  equations  in  the  two  unknowns,  w°  and 
Zi.  At  this  point  the  complementary  beam  ( including  the  mode  coeffi- 
cient, y ) is  completely  specified.  Within  the  transmitter  plane  the 
total  field  is 


E(x,y,0)  = Eg(x,y,0) 


(F— 14) 
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APPENDIX  G 


COMPUTER  PROGRAM  LISTINGS 


This  appendix  contains  listings  and  brief  discussions  of  the  com- 
puter programs  used  to  evaluate  the  temporal  scintillation  spectrum  for 
the  clean  gauss Ian  beam  (Chapter  III)  and  the  dirty  gauss Ian  beam 
(Chapter  IV)  for  homogeneous  and/or  localized  turbulence.  These  pro- 
grams are  designed  to  evaluate  the  expression 


SDU) 


[ dz  [ die  f (z)[w/w0)2  + 
J0  -00 


where 


i2(o,o,l)c2  * i2(o,o,l)  Jo  J — 

d * 

+ (1.077/LJ)2r11/6  x |*(z,  j-*  <.w)  - H*(z.-  (G-la) 

hL  S-  " At) 

V uo  J 


x Ex 


(x, 


v^XL”(  1 -Z  )oj/u)( 

ZT~ 


’0  /xF(1-z)k\ 

- — n ) E*<xo-V  • 


(G-lb) 


C is  a constant  defined  as 


c = 256/2r 


^ (-cosjj-.)  ■ 48.527 


(6-lc) 


and  f(z)  Is  a longitudinal  turbulence  strength  distribution  function. 
For  homogeneous  turbulence,  as  In  program  one. 


f(z)  - 1 . (G-2) 

In  program  two  we  used 

f (z)  * 100  e"100z  . (G-3) 
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Equation  (6-lb)  mbs  obtained  from  Equation  (89b)  In  the  text  by  assuming 
a point  detector i 


w(r)  ■ «(|r-r0|)  . . 

For  the  clean  gausslan  beam  the  field,  E,  In  Equation  (6-lb)  was  given 
by  Equation  (194)  In  the  text.  The  field  of  the  dirty  gausslan  beam 
was  given  by  Equation  (202b). 


To  simplify  the  programming,  use  was  made  of  the  change  of  variable 

' K 


The  numerical  algorithm  employed  in  approximating  the  integrals 
in  Equation  (G-4)  was  Gaussian  quadrature  [95,96,97].  This  technique 
basically  involves  fitting  a Hermite  polynomial  to  the  function  to  be 
integrated  and  integrating  the  polynomial  over  the  desired  limits.  To 
approximate  the  integration  over  the  range  variable  z,  a single  32 
point  Gaussian  quadrature  routine  was  used.  The  32  point  routine  as 
opposed  to  a higher  order  quadrature  formula  was  deemed  adequate  be- 
cause the  z integrand  was  slowly  varying.  The  k integration  was  approx- 
imated by  breaking  the  range  into  segments  thus 


co  M 

d<  F(k)  * I 
0 i=l 


i.A 

die  F(k) 

( i -1 ) A 


(6-5) 


and  employing  an  eight  point  Gaussian  quadrature  routine  on  each  seg- 
ment. In  Equation  (G-5)  the  width,  a,  of  the  segments  was  chosen  such 
that 


(6-6) 


and  the  upper  limit,  M,  on  the  sum  was  such  that  the  last  segment 
integrated  made  a negligible  contribution  to  the  integral,  i.e., 

(Ma 

dtcF(ic) 

J(M-1)A 

M-l  fiA 
I d<F( 

i=l  ‘ ( i -1 ) A 


This  method  of  breaking  a large  integration  range  into  smaller  pieces 
is  commonly  called  composite  or  compound  integration  [98].  The  advan- 
tages of  this  technique  are  a smaller  error  and  the  ability  to  approxi- 
mate integrations  over  infinite  or  semi -infinite  ranges  [98]. 

Listings  of  the  programs  used  in  Chapters  III  and  IV  follow.  Also 
included  are  listings  for  the  eight  and  32  point  Gaussian  quadrature 
subroutines.  The  notation  used  in  programs  one  and  two  is  consistent 
with  that  of  the  text  and  of  this  appendix. 
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Program  one  calculates  the  dirty  beam  scintillation  spectrum  for 
homogeneous  turbulence.  To  obtain  the  scintillation  spectrum  of  a 
clean  beam  from  this  program,  the  complex  variable  C2  is  set  to  zero. 
The  dirty  beam  spectrum  for  localized  turbulence  is  calculated  by 
program  two.  A linear  combination  of  the  results  of  these  two  pro- 
grams yields  the  dirty  beam  spectrum  for  the  combination  of  homo- 
geneous and  localized  turbulence. 


PROGRAM  ONE 


BEST  AVAILABLE  COPY 


C SCINTILLATION  SPECTRUM  FOR  DIRTY  BEAM  (PHASE  OBJECT) 

C WITH  OFF-AXIS  POINT  DETECTOR  AND  LOCALIZED  TURBULENCE. 

C ZERO-ZERO  ORDER  APPROXIMATION  TO  DIRT  SPOT 
OPTIONS  DP 
INCLUDE  DQG88»DOH83 
EXTERNAL  ZINT.EFI ELD 
REAL  IOSO 

E COMPLEX  A,B, 8C, AL.BL.CI ,C2, EFO.EFIELD 

COMMON  AL.3L.CI ,C2,OMEG,FOMEG I .FOMEG2 ,F 1 Z.F2Z, EFO, EKSO, WYEO 
I.PI.FRES 

CALL  A SSI GNC6HTMPOUT, O. ,6 ) 

CALL  FERR(O) 

CALL  DEASSH 
PI =4 . *ATAN ( I . ) 

C POSITION  OF  POINT  DETECTOR 
EKS0=0. 

WYE0=SQRT(.2> 

C DIAMETER  AND  PHASE  OF  DIRT  SPOT 
DEL-.9I7E-J 
PHI »5.333580778 

C WAIST  SPOTSIZES  AND  LOCATIONS  FOR- MAIN  AND  COMPLEMENTARY  BEAMS 
W0-20.2E-6 
Wl»3. I46345943E-4 
ZO-I  . 

Z I “.4080236 1 59 

C RANGE,  OUTER  SCALE,  AND  WAVELENGTH 
EL-I.E3 
ELO“l . 

ALAM“ . 6328E-6 
FRES=SORT  ( ALA.M*EL ) 

FELOP“ I .077/(EL0/FRES) 

C AA-N/L 

AA= ALAM/C  PI *WO*WO) 

C PB-M/L 

BB= ALAM/C  PI *W 1 *WI  ) 

C A“AL?HA(ZO> 

A“WO*WO*CMPLX( I . , AA*ZO) 

A=I ./A 
C B=BETA(ZI ) 

B=W  I *WI *CMPLX  ( I . , EB*Z  I ) 

B=l ./B 
BC*CONJG(B) 

C AL=ALPHA( L+ZO ) 

AL=WO*WO*CMPLX< I . , AA*(ZO+EL > ) 

AL-I ./AL 

C Cl “ALPHA (L+0)*WO 
C1=WC*AL 

C I 0SC= I (0, 0, L) +*2 

I OSO=REAL(C I *CONJG( Cl ))**2 
C BL”8ETA( L+Z I ) 

8L-W I *rt I *CMPLX ( I. ,BB*(ZI+EL)) 

BL»I ./3L 

C C2-EETA(L+ZI)*GAAMA00*WI 

C2“( I .-CEXPC-C  A+BC ) *( DEL/2. )**2 ) )*( I .-CEXP(CMPLX(0. , PHI J J > 
C2»C2*A*3C/( A+BC) 

C2-2 . *W I *W0*C2 
C2“C2*BL*WI 
EFO“EFIELD( EKSO, WYEO) 

C OMEG-OMEGA/OMEGAO 
OMEG-I .E-5 
FRIMC«IO.**< | ,/8. ) 

CONST-48.52728232 
DO  50  1-1,57 
FOMEG I « (OMEG) **2 

FOMPOjaPoupni np**5 
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BEST  AVAILABLE  COPY 


■ calit  DdC32(  o. ; rv;zrNT , ans  ) 

ANS-. 5*ANS*CONST/I0SQ 
WRITECo,- )OMEG, ANS 
50  OMEG-OMEG*FRINC 
CALL  EXIT 
END 
C 

C THIS  IS  THE  INTEGRAND  FOR  THE  PATH  INTEGRATION 
FUNCTION  ZINT(Z) 

EXTERNAL  AKINT 
COMPLEX  AL, BL.CI .C2.EF0 

COMMON  AL,BL,CI,C2,OMEG,FOMEGI,FOMEC2.FIZ.F2Z.EFO.EXSO,WYEO 
I.PI.FRES 

FIZ-(I.-Z)/(2.*PI) 

F2Z-FIZ*FRES 

AUPR-F0MEG2 

IFCFOMEGI .GE. I . JAUPR- I ./SCRTTFOMEGI ) 

Z I NT-0. 

ALWR-O. 

5 CALL  DQG8(ALhR.AUPR, AKINT. ANS) 

ZINT-ZINT+ANS 

ALNR-AUPR 

AUPR-2.*AUPR 

I FC  ANS. LT.  I . E-5*Z  IN'T)  RETURN 
C0T05 
END 
C 

C THIS  IS  THE  KAPPA  INTEGRAND 
FUNCTION  AKINT ( AKP) 

COMPLEX  AL, BL.CI , C2 , EFO.FCN, EFI ELD, F6.F7 

COMMON  AL, BL.CI ,C2,0MEG,F0MEG1 .FOME02 ,F I Z ,F2Z,EF0, EKSO, WYEO 

F l-AKP+FOMEGl 

F2-AKP+FOMEG2 

F3-FI Z*F1 *.5 

F 4»F2Z*S0RT ( AKP ) 

F5«F2Z*OMEG 

F<5«EFO*CEXP(CMPLX(O..F3)) 

F7-CONJG(F6 ) 

FCN-EFI ELD( EKS0-F5, WYEO-F  4 ) *F7 
FCN-FCN-CONJG ( EF I ELD ( EK  S0*F5 , WYEO+F  4 ) >* F6 
AKINT-REAL(FCN*CONJG( FCN ) ) 

FCN-EF I ELD  ( EKS0-F5 , WYE0+F4 ) *F7 
FCN-FCN-CONJG (EFI ELD( EKS0+F5, WYE0-F4 ) ) * F6 
AKI NT-AKI NT+REAL( FCN*a)NJG( FCN) ) 

AKI NT- AKI NT*F2 **( -I .8  3333333333 ) * AKP**< - . 5) 

RETURN 

END 

C 

C THIS  FUNCTION  CALCULATES  THE  EFIELD  AT  THE  RECEIVER  PLANE. 
COMPLEX  FUNCTION  EFIELD(X.Y) 

COMPLEX  AL, BL.CI ,C2 
COMMON  AL. BL.CI ,C2 
RSQ-X*X*Y*Y 

EFI ELD-C I *CEXP ( -AL*RSO ) -C2*CEXP < -BL*RSQ ) 

RETURN 


BESI'AVAItABLE  COPY 


PROGRAM  TWO 


C SCINTILLATION  SPECTRUM  FOR  DIRTY  BEAM  (PHASE  OBJECT) 

C WITH  OFF-AXIS  POINT  DETECTOR. 

C ZERO-ZERO  ORDER  APPROXIMATION  TO  DIRT  SPOT 
OPTIONS  DP 
INCLUDE  DQG8E IDQG32B 
EXTERNAL  ZINT.EFIELD 
REAL  I OSQ 

COMPLEX  A,3,BC,AL,BL.CI ,C2,EFO,EFIELD 

COMMON  AL , BL,  C I , C2 , OMEG . FO.MEG  I , FOMEG2 , F 1 Z , F2Z,  EFO,  EK SO,  WYEO 
I.PI.FRES 

CALL  ASSIGN(6HTMP0UT,0. .6) 

CALL  DEASSN 
CALL  FERR(O) 

PI-4.*ATAN(I .) 

C POSITION  OF  POINT  DETECTOR 
EKSO-O. 

WYEO=SQRT( .4 ) 

C DIAMETER  AND  PHASE  OF  DIRT  SPOT 
Da-.9f  7E-3 
PHI -5. 333580773 

C WAIST  SPOTSIZES  AND  LOCATIONS  FOR' MAIN  AND  COMPLEMENTARY  BEAMS 
W0-20.2E-6 
hi  m3. I46345943E-4 
ZO-I. 

Z I *.4080236 1 59 

C RANGE,  OUTER  SCALE,  AND  WAVELENGTH 
EL-I.E3 
ELO- ! .E2 
ALA.M-.6328E-6 
FRES-SQRT  ( A.LAM*EL ) 

FBLOPm 1.077/1 ELO/FRES ) 

C AA-N/L 

AA-ALAM/(PI*WO*WO) 

C 93-M/L 

BB«ALAM/(PI*WI*WI ) 

C A-ALPHA(ZO) 

A-HO*WO*C.MPLX(l  . , AA*ZO) 

A-l ./A 
C B-BETA(ZO) 

3 »W I *W I *CMPLX  (1  . , BB*Z I ) 

3- I ./B 
BC-CONJG(B) 

C AL»ALPHA(L+ZO) 

AL-WC*WO*CMPLX( I . , AA*<ZO+EL ) ) 

AL-I ./AL 

C C1-ALPHA(L+Z0)*W0 
C I -WO*AL 

C I OSC-I ( 0,0,L) **2 

I0SQ*REAL<CI*C0NJG(C1 ))**2 
C BL-EETACL+ZI ) 

BL-WI*WI*CMPLX< I .,BB*(ZI+EL)> 

3L-I./BL 

C C2-BETA(L+ZI )*GAMMAOO*WI 

C2-(l .-CEXP(-(A+BC)*(DEL/2.)**2))*( l.-CEXP(CMPLXCO. ,PHI ))) 
C2-C2*A*BC/ ( A+BC ) 

C2-2.*WI*W0*C2 

C2-C2*BL*WI 

EFO«EFIELD(EKSO,WYEO) 

C OMEC-OMSGA/OMEGAO 
OMEG- I . E-5 
FRINC-I 0. **( 1 ,/8. ) 

CONST -48. 52723232 
DO  50  1-1,57 
FOMEGI -( OMEG) **2 
FOMEG2-FOMEGI+FELOP—2 
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BESt  AVAILABLE  COPY 


CALL  D0H8(0.,.2,ZINT.ANS) 

A NS- . 5*ANS*CONST/ I OSO 
ANS-ANS*EL/IO. 

WR I TE ( 6 )0  MEG , ANS 
50  OMEG«OMEG*FRINC 
CALL  EXIT 
END 
C 

C THIS  IS  THE  INTEGRAND  FOR  THE  PATH  INTEGRATION 
FUNCTION  ZINT(Z) 

EXTERNAL  AKINT 
COMPLEX  AL.BL.C  I,  C2.EF0 

COMMON  AL.BL.C) .C2.0MEG, F0MEG1 , F0MEG2 ,F IZ .F2Z .EFO.EKSO, WYEO 
1.PI.FRE5 

F1Z*(I.-Z)/C2.*PI) 

F2Z-FIZ*FRES 

AUPR"F0MEG2 

IF(FOMEG) .GE. I. > AUPR=  I ./SQRTC  FOMEGI ) 

ZINT-O. 

ALWR=0. 

5 CALL  DQG8<  ALP-R.AUPR, AKINT, ANS) 

ZINT=ZINT+ANS 

ALWR«AUPR 

AUPR=2.*AUPR 

I F ( ANS. LT . 1 . E-5*Z I NT) GOTO  1 0 
GOTO  5 

10  Z INT=ZINT*EXP(-Z* 1 00. ) 

RETURN 

END 

C 

C THIS  IS  THE  KAPPA  INTEGRAND 
FUNCTION  AKINT(AKP) 

COMPLEX  AL.BL.C) ,C2 .EFO.FCN ,EFI ELD.F6 ,F7 

COMMON  AL, BL,C I ,C2 .OMEG.FOMEG I .F0MEG2 ,F IZ.F2Z.EF0, EKSO, WYEO 

F1»AKP+F0MEGI 

F2=AKP+FOMEG2 

F3=F 1 Z*F 1 *.5 

F/=F2Z*S0RT(AKP) 

F5«F2Z*OMEG 

F6=EF0*CEXP(CMPLX (0. , F3 ) ) 

F7*C0NJG(F6) 

FCN»EFIELD(EKS0-F5,WYE0-F4)*F7 
FCN=FCN-CON  JG  < EF  I ELD  C EK  SCH-F  5 . NYEO+F  4 ) >*F6 
AKI NT-REALC FCN*CONJG( FCN ) ) 

FCN-  EF  I ELD( EKS0-F5 . KYE0+F4 ) *F7 
FCN-FCN-CON JG ( EF I ELD ( EKS0+F5 . WYE0-F4 ) )*F6 
AKINT-AKINT+REAL(FCN*CONJG(FCN) ) 

AKINT«AKINT*F2**<-1 .83333333333 )*AKP**(-. 5) 

RETURN 

END 

C 

C THIS  FUNCTION  CALCULATES  THE  EFIELD  AT  THE  RECEIVER  PLANE. 
COMPLEX  FUNCTION  EFIELD(X.Y) 

COMPLEX  AL.BL.CI ,C2 
COMMON  AL.BL.C) ,C2 
RSO-X*X+Y*Y 

EFIELD-CI*CEXPC-AL*RS0)-C2*CEXP<-BL*RS0) 

RETURN 

END 


nn  ooooooooooooooooooooooooooooooooooo 


BEST  AVAILABLE  COPY 

EIGHT  POINT  GAUSSIAN  QUADRATURE  SUBROUTINE 


12 

SUBROUTINE  DQG8  13 

14 

PURPOSE  15 

TO  COMPUTE  INTEGRALY FCTY  X) , SUMMED  OVER  X FROM  XL  TO  XU)  16 

17 

USAGE  18 

CALL  DQG8  ( XL ,XU, FCT, Y)  19 

PARAMETER  FCT  REQUIRES  AN  EXTERNAL  STATEMENT  20 

21 

DESCRIPTION  OF  PARAMETERS  22 

XL  -DOUBLE  PRECISION  LOWER  BOUND  OF  THE  INTERVAL.  23 

XU  -DOUBLE  PRECISION  UPPER  BOUND  OF  THE  INTERVAL.  24 

FCT  -THE  NAME  OF  AN  EXTERNAL  DOUBLE  PRECISION  FUNCTION  25 

SUBPROGRAM  USED.  26 

Y -THE  RESULTING  DOUBLE  PRECISION  INTEGRAL  VALUE.  27 

28 

REMARKS  29 

NONE  30 

31 

SUBROUTINES  AND  FUNCTION  SUBPROGRAMS  REQUIRED  32 

THE  EXTERNAL  DOUBLE  PRECISION  FUNCTION  SUBPROGRAM  FCT(X)  33 
MUST  BE  FURNISHED  BY  THE  USER.  34 

35 

METHOD  36 

EVALUATION  IS  DONE  BY  MEANS  OF  8-POINT  GAUSS  QUADRATURE  37 

FORMULA,  WHICH  INTEGRATES  POLYNOMIALS  UP  TO  DEGREE  15  38 

EXACTLY#  FOR  REFERENCE  SEE  39 

HANDBOOK  OF  MATHEMATICAL  FUNCTIONS  WITH  FORMULAS,  GRAPHS.  AND  40 

MATHEMATICAL  TABLES, EDITED  BY  M.  ABRAMOWITZ  AND  I.  A.  STEGUN 
U.S.  GOVERNMENT  PRINTING  OFFICE,  WASHINGTON  D.C.  42 

PAGES  887,888,919  43 

44 

46 

47 

SUBROUTINE  DQG8  ( XL, XU, FCT, Y)  48 

49 

50 

DOUBLE  PRECISION  XL,  XU,  Y,  A,  B,  C,  FCT  51 

C 52 

A«. EDO* (XU* XL)  53 

B-XU-XL  54 

C-.480I4492824876D0+B  55 

Y=. 506 1 42  68 1 45 1 88D- 1 * (FCT ( A+C )+FCT( A-C ) ) 56 

C*. 3983332387068 1 *8  57 

Y-Y+. I ! 1 1 905I722668D0*(FCT( A+C)+FCT( A-C  ) ) 58 

C». 262766204958 I6*B  59 

Y-Y*.  I5685332293894D0+Y FCT( A*C)*FCT( A-C ) ) 60 

C*.91 7 1 732 1 24782D-I *3  61 

Y«B*(Y+. 181 34 1891 6891 8D0* (FCTY A+C >+FCT( A-C)))  62 

RETURN  63 

END  64 

• i 
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1 


best  available  copy 


THIRTY-TWO  POINT  GAUSSIAN  QUADRATURE  SUBROUTINE 


SUBROUTINE  DCG32 


PURPOSE 

TO  COMPUTE  INTEGRAL!  FCT(  X) , SUMMED  OVER  X FROM  XL  TO  X'J) 


USAGE 

CALL  DGG32  ( XL, XU, FCT, Y) 

PARAMETER  FCT  REQUIRES  AN  EXTERNAL  STATEMENT 


DESCRIPTION  OF  PARAMETERS 


XL 

XU 

FCT 


-DOUBLE  PRECISION  LOWER  BOUND  OF  THE  INTERVAL. 
-DOUBLE  PRECISION  UPPER  BOUND  OF  THE  INTERVAL. 
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15 
10 

17 

18 

19 

20 
21 
22 

23 

24 


-THE  NAME  OF  AN  EXTERNAL  DOUBLE  PRECISION  FUNCTION  25 


SUBPROGRAM  USED. 

-THE  RESULTING  DOUBLE  PRECISION  INTEGRAL  VALUE. 


REMARKS 

NONE 


SUBROUTINES  AND  FUNCTION  SUBPROGRAMS  REQUIRED 

THE  EXTERNAL  DOUBLE  PRECISION  FUNCTION  SUBPROGRAM  FCT(X) 
MUST  BE  FURNISHED  BY  THE  USER. 


METHOD 

EVALUATION  IS  DONE  BY  MEANS  OF  32-POINT  GAUSS  QUADRATURE 
FORMULA,  WHICH  INTEGRATES  POLYNOMIALS  UP  TO  DEGREE  63 
EXACTLY.  FOR  REFERENCE,  SEE 
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39 


HANDBOOK  OF  MATHEMATICAL  FUNCTIONS  WITH  FORMULAS,  GRAPHS,  AND40 


MATHEMATICAL  TABLES, EDITED  BY  M.  ABR AM0WIT2  AND  I.  A.  STEGUN 
U.S.  GOVERNMENT  PRINTING  OFFICE,  WASHINGTON  D.C. 

PAGES  887, 888,919 


SUBROUTINE  D0G32C  XL,XU,FCT, Y) 


DOUBLE  PRECISION  XL,  XU.  Y,  A,  B,  C,  FCT 


A-.5D0*(XU+XL> 

B-XU-XL 

C*. 49863 1 93092474C0*B 

Y» . 3509305004735 OD-2*  C FCT  C A+C ) +FCTC  A -C  > ) 
C*. 492805 75577263D0*B 

Y-Y+. 8 1 37 1 97365452D-2*( FCTC  A+C) +FCTC  A-C ) ) 
C=. 482381 I2779375D0*B 

Y-Y+. 1269603265463! D- l*(FCT( A+C  J+FCTC  A-C) ) 
C«.46745303796886D0*B 

Y «Y+ . 1 7 1 3 693 1 4 5 65 1 0D- 1 * ( FCT ( A+C  > + FCT ( A-C  > ) 
C-.448l60577e8302D0*B 

Y»Y+.2 1 4 1 79490 1 I I 1 3D- 1* (FCTC  A+C ) ♦FCTC A-C ) ) 
C*.42468380686628DO*B 

Y-Y-t-. 25499029631 I 88D- I* (FCTC A+C J+FCTC  A-C) ) 
C-.39724 I 89798397D0*B 

Y-Y+.29342046739267D- 1* (FCTC A+C J+FCTC A-C ) ) 
C-. 36609 1 059270 1 4D0*B 

Y-Y+.329I 1 1 I 13881 80D- 1* (FCTC  A+C J+FCTC A-C) ) 
C*.  33 1 522 1 33465 1 0D0*B 

Y»Y+ . 36 1 72897054424D- I*  < FCTC  A+C ) +FCTC A-C ) ) 
C«.29385787862038D0*B 

Y-Y+ . 390969478935  35D- 1 * ( FCTC  A+C ) *FCT ( A-C ) ) 
C-.  25344995 4466 1 1 D0*B 

Y-Y+ . 4 1 655962  I 1 3473 D- 1 * C FCTC  A*C J+FCTC  A-C ) ) 
C-. 2 1 06 7563806531  D0*B 

Y-Y+ . 438260465022  C I D- 1 * ( FCTC  A+C ) ♦FCTC A-C ) ) 
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RBI  "AVAILABLE  COPY 


C-.  165934301  I 4 1 06DO*B 

Y“Y*. 4558693534788 I D- I *( FCT( A+C )+FCT( A-C ) ) 

O.  I I 964368  1 1 2606 D0*B 

Y-Y+-.  46922 1 955404 C2D-  l*(FCT(A*C)*FCT(  A-C ) ) 

O. 72235980791 39D-I*B 
Y«Y«-.478I936C039637D-I*(FCT(A*C)+FCT(  A-C)) 

C». 24 1 53832843869D- 1*9 

Y»B*(  Y+. 4827004425736 3D- 1 *(  FCT(  A+O+FCTC  A-C  ) ) ) 

RETURN 

END 
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